{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Pima Indians Diabetes Data Set\n",
    "数据说明： Pima Indians Diabetes Data Set（皮马印第安人糖尿病数据集） 根据现有的医疗信息预测5年内皮马印第安人糖尿病发作的概率。\n",
    "\n",
    "数据集共9个字段: <br>\n",
    "第0列为怀孕次数； <br>\n",
    "第1列为口服葡萄糖耐量试验中2小时后的血浆葡萄糖浓度；<br> \n",
    "第2列为舒张压（单位:mm Hg） <br>\n",
    "第3列为三头肌皮褶厚度（单位：mm） <br>\n",
    "第4列为餐后血清胰岛素（单位:mm） <br>\n",
    "第5列为体重指数（体重（公斤）/ 身高（米）*2） <br>\n",
    "第6列为糖尿病家系作用 <br>\n",
    "第7列为年龄 <br>\n",
    "第8列为分类变量（0或1）<br>\n",
    "\n",
    "导入标准化后的数据，进行logicRegression回归训练模型，分别用不同的模型进行测试并验证log loss等"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "#import必要的工具包，用于文件读取／特征编码\n",
    "\n",
    "import pandas as pd\n",
    "import numpy as np\n",
    "from sklearn.model_selection import GridSearchCV\n",
    "\n",
    "import matplotlib.pyplot as plt\n",
    "%matplotlib inline\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>pregnants</th>\n",
       "      <th>Plasma_glucose_concentration</th>\n",
       "      <th>blood_pressure</th>\n",
       "      <th>Triceps_skin_fold_thickness</th>\n",
       "      <th>serum_insulin</th>\n",
       "      <th>BMI</th>\n",
       "      <th>Diabetes_pedigree_function</th>\n",
       "      <th>Age</th>\n",
       "      <th>Triceps_skin_fold_thickness_Missing</th>\n",
       "      <th>Target</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>0.639947</td>\n",
       "      <td>0.866045</td>\n",
       "      <td>-0.031990</td>\n",
       "      <td>0.670643</td>\n",
       "      <td>-0.181541</td>\n",
       "      <td>0.166619</td>\n",
       "      <td>0.468492</td>\n",
       "      <td>1.425995</td>\n",
       "      <td>-0.647760</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>-0.844885</td>\n",
       "      <td>-1.205066</td>\n",
       "      <td>-0.528319</td>\n",
       "      <td>-0.012301</td>\n",
       "      <td>-0.181541</td>\n",
       "      <td>-0.852200</td>\n",
       "      <td>-0.365061</td>\n",
       "      <td>-0.190672</td>\n",
       "      <td>-0.647760</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>1.233880</td>\n",
       "      <td>2.016662</td>\n",
       "      <td>-0.693761</td>\n",
       "      <td>-0.012301</td>\n",
       "      <td>-0.181541</td>\n",
       "      <td>-1.332500</td>\n",
       "      <td>0.604397</td>\n",
       "      <td>-0.105584</td>\n",
       "      <td>1.543781</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>-0.844885</td>\n",
       "      <td>-1.073567</td>\n",
       "      <td>-0.528319</td>\n",
       "      <td>-0.695245</td>\n",
       "      <td>-0.540642</td>\n",
       "      <td>-0.633881</td>\n",
       "      <td>-0.920763</td>\n",
       "      <td>-1.041549</td>\n",
       "      <td>-0.647760</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>-1.141852</td>\n",
       "      <td>0.504422</td>\n",
       "      <td>-2.679076</td>\n",
       "      <td>0.670643</td>\n",
       "      <td>0.316566</td>\n",
       "      <td>1.549303</td>\n",
       "      <td>5.484909</td>\n",
       "      <td>-0.020496</td>\n",
       "      <td>-0.647760</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "   pregnants  Plasma_glucose_concentration  blood_pressure  \\\n",
       "0   0.639947                      0.866045       -0.031990   \n",
       "1  -0.844885                     -1.205066       -0.528319   \n",
       "2   1.233880                      2.016662       -0.693761   \n",
       "3  -0.844885                     -1.073567       -0.528319   \n",
       "4  -1.141852                      0.504422       -2.679076   \n",
       "\n",
       "   Triceps_skin_fold_thickness  serum_insulin       BMI  \\\n",
       "0                     0.670643      -0.181541  0.166619   \n",
       "1                    -0.012301      -0.181541 -0.852200   \n",
       "2                    -0.012301      -0.181541 -1.332500   \n",
       "3                    -0.695245      -0.540642 -0.633881   \n",
       "4                     0.670643       0.316566  1.549303   \n",
       "\n",
       "   Diabetes_pedigree_function       Age  Triceps_skin_fold_thickness_Missing  \\\n",
       "0                    0.468492  1.425995                            -0.647760   \n",
       "1                   -0.365061 -0.190672                            -0.647760   \n",
       "2                    0.604397 -0.105584                             1.543781   \n",
       "3                   -0.920763 -1.041549                            -0.647760   \n",
       "4                    5.484909 -0.020496                            -0.647760   \n",
       "\n",
       "   Target  \n",
       "0       1  \n",
       "1       0  \n",
       "2       1  \n",
       "3       0  \n",
       "4       1  "
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "#使用数据为在前面处理做特征标准化整理后的数据\n",
    "train = pd.read_csv(\"FE_pima-indians-diabetes.csv\")\n",
    "train.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "y_train = train['Target']   \n",
    "X_train = train.drop([\"Target\"], axis=1)\n",
    "\n",
    "#保存特征名字以备后用（可视化）\n",
    "feat_names = X_train.columns "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 默认参数的Logistic Regression"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.linear_model import LogisticRegression\n",
    "lr = LogisticRegression()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "logloss of each fold is:  [0.48759712 0.52695217 0.45678156 0.42161373 0.48250956]\n",
      "cv logloss is: 0.47509082884510445\n",
      "accuracy of each fold is:  [0.75974026 0.74025974 0.76623377 0.79738562 0.76470588]\n",
      "cv accuracy is: 0.7656650539003481\n"
     ]
    }
   ],
   "source": [
    "# import warnings filter\n",
    "from warnings import simplefilter\n",
    "# ignore all future warnings\n",
    "simplefilter(action='ignore', category=FutureWarning)\n",
    "\n",
    "# 交叉验证用于评估模型性能和进行参数调优（模型选择）\n",
    "#分类任务中交叉验证缺省是采用StratifiedKFold\n",
    "#数据集比较大，采用5折交叉验证\n",
    "from sklearn.model_selection import cross_val_score\n",
    "loss = cross_val_score(lr, X_train, y_train, cv=5, scoring='neg_log_loss')\n",
    "#%timeit loss_sparse = cross_val_score(lr, X_train_sparse, y_train, cv=3, scoring='neg_log_loss')\n",
    "print('logloss of each fold is: ',-loss)\n",
    "print('cv logloss is:', -loss.mean())\n",
    "\n",
    "accuracy = cross_val_score(lr, X_train, y_train, cv=5, scoring='accuracy')\n",
    "print('accuracy of each fold is: ',accuracy)\n",
    "print('cv accuracy is:', accuracy.mean())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Logistic Regression + GridSearchCV\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "#为了比较GridSearchCV和LogisticRegressionCV，两者用相同的交叉验证数据分割\n",
    "from sklearn.model_selection import StratifiedKFold\n",
    "fold = StratifiedKFold(n_splits=5, shuffle=True, random_state=777)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "GridSearchCV(cv=StratifiedKFold(n_splits=5, random_state=777, shuffle=True),\n",
       "       error_score='raise-deprecating',\n",
       "       estimator=LogisticRegression(C=1.0, class_weight=None, dual=False, fit_intercept=True,\n",
       "          intercept_scaling=1, max_iter=100, multi_class='warn',\n",
       "          n_jobs=None, penalty='l2', random_state=None, solver='warn',\n",
       "          tol=0.0001, verbose=0, warm_start=False),\n",
       "       fit_params=None, iid='warn', n_jobs=None,\n",
       "       param_grid={'penalty': ['l1', 'l2'], 'C': [0.001, 0.01, 0.1, 1, 10, 100, 1000]},\n",
       "       pre_dispatch='2*n_jobs', refit=True, return_train_score='warn',\n",
       "       scoring='neg_log_loss', verbose=0)"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "from sklearn.model_selection import GridSearchCV\n",
    "from sklearn.linear_model import LogisticRegression\n",
    "\n",
    "#需要调优的参数\n",
    "# 请尝试将L1正则和L2正则分开，并配合合适的优化求解算法（slover）\n",
    "#tuned_parameters = {'penalty':['l1','l2'],\n",
    "#                   'C': [0.001, 0.01, 0.1, 1, 10, 100, 1000]\n",
    "#                   }\n",
    "penaltys = ['l1','l2']\n",
    "Cs = [0.001, 0.01, 0.1, 1, 10, 100, 1000]\n",
    "tuned_parameters = dict(penalty = penaltys, C = Cs)\n",
    "#LogisticRegression 的 penalty默认为L2\n",
    "lr_penalty= LogisticRegression()\n",
    "grid= GridSearchCV(lr_penalty, tuned_parameters,cv=fold, scoring='neg_log_loss')\n",
    "grid.fit(X_train,y_train)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.47529906721481324\n",
      "{'C': 1, 'penalty': 'l1'}\n"
     ]
    }
   ],
   "source": [
    "# examine the best model\n",
    "print(-grid.best_score_)\n",
    "print(grid.best_params_)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYsAAAEKCAYAAADjDHn2AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzt3Xd8VFX+//HXZ2bSSUggAWkSQFBaQgmgYoNVREB6Fym7iK403dVd1++udV3bT0UURUFZZVFEuoiAiIiilNB7bwFCCTWkZ87vjxkwQCCTMrmZ5PPcxzwy9845M+8b3Pnk3nPvuWKMQSmllLoem9UBlFJKlXxaLJRSSuVJi4VSSqk8abFQSimVJy0WSiml8qTFQimlVJ60WCillMqTFgullFJ50mKhlFIqTw6rAxSVyMhIEx0dbXUMpZTyKWvWrDlpjInKq51Xi4WItAfeBezARGPMa1e8/g7Qxr0YDFQyxoS7XxsE/NP92r+NMZ9d77Oio6OJj48vyvhKKVXqicgBT9p5rViIiB0YB9wHJACrRWSuMWbrxTbGmCdztB8JNHU/rwA8D8QBBljj7nvaW3mVUkpdmzfHLFoCu40xe40xGcBUoMt12vcDvnQ/vx/43hhzyl0gvgfaezGrUkqp6/BmsagGHMqxnOBedxURqQnUApbkp6+IDBOReBGJP3HiRJGEVkopdbWSMsDdF5hujMnOTydjzMfAxwBxcXE617pSZVRmZiYJCQmkpaVZHaXECgwMpHr16vj5+RWovzeLxWGgRo7l6u51uekLDL+i7z1X9F1ahNmUUqVIQkICoaGhREdHIyJWxylxjDEkJSWRkJBArVq1CvQe3jwMtRqoKyK1RMQfV0GYe2UjEbkFiAB+y7F6IdBORCJEJAJo516nlFJXSUtLo2LFiloorkFEqFixYqH2vLy2Z2GMyRKREbi+5O3Ap8aYLSLyEhBvjLlYOPoCU02OW/YZY06JyMu4Cg7AS8aYU97KqpTyfVoorq+wvx+vjlkYY+YD869Y99wVyy9co++nwKdeC+eWnZXF6omjqPnAaKrUvNnbH6eUKiH6fOQ6mPHVo7dZnMQ3lPnpPo7s30aDxNnIpI4k7NlidRyllI8qV67cpeft27cnPDycTp065dp2+PDhNGnShAYNGhAUFESTJk1o0qQJ06dPz9dnrl27lgULFhQqt6fKfLGocVNjTnT/mgDS8J/cif3b11kdSSnl455++mkmT558zdfHjRvH+vXrmT9/PnXq1GH9+vWsX7+enj175utztFgUszoxrTnXZzZ2nIRO7cLuzSutjqSU8mF/+MMfCA0NLVDfXbt2cf/999O8eXPuuusudu7cCcDUqVNp1KgRsbGxtGnThtTUVF566SWmTJlSoL2S/Cop11lYrmb9OA499A0BU7oQOb072zO+5JZmd1kdSymVTy9+s4WtR87l2W7rUVebi2MX19OgahjPP9iw0Nk8MWzYMCZOnEidOnVYvnw5I0aMYNGiRbz44ossXbqUypUrc+bMGYKCgnjuuefYvHkzY8aM8XouLRY51Kgbw9HB80n9rDNV5/RhS+bnNGx1n9WxlFJlxJkzZ1ixYgU9evS4tC4rKwuA1q1bM3DgQHr16kX37t2LPZsWiytUia7PiaELOPdJR6LnD2BDxifE3pn7IJVSquTxdA+gJJ4NZYwhMjKS9evXX/XahAkTWLlyJfPmzaNZs2asW1e846s6ZpGLqGp1CB62iJP2StRbPIQ1S2ZYHUkpVQZERERQpUoVZs2aBYDT6WTDhg0A7N27l1tvvZWXX36ZiIgIDh8+TGhoKOfPny+WbFosrqHCDTcS8fgijjqq0/inYaxcMMXqSEopH3HnnXfSq1cvfvjhB6pXr87ChZ5PQDF16lTGjx9PbGwsDRs2ZN68eQA8+eSTNG7cmMaNG9OmTRsaNWpE27Zt2bBhA02bNtUBbiuFRVZBRizk4AcdafbbSH7LTOO2B/9kdSylVAmUnJx86fnPP//sUZ/o6Gg2b9582bratWvnWlzmzr1qtiSioqKK7aZvWizyEBpRCcfIRex7vwMt4//KL5np3NH9catjKaUKqSSNVfgCPQzlgaCwCG4cvYCdQbHcvuFZlk19y+pISilVrLRYeCgwpDw3PTGfbSEtuGv7Syyd/Ao55j5USqlSTYtFPvgFhnDz6LlsCr2De/a8wY+TntOCoZQqE7RY5JMjIIiGo2axsXxb2h4cy48fP43TqQVDKVW6abEoAJufP41HTWNjxQdoe3QCSz8cSXa20+pYSqn8mNTR9VAe0WJRQGL3o/HwKWys3I22Jyaz7P1hZGXl6xbiSqlSpLinKJ81axZvvvlmoXN7Sk+dLQSx2Yl5bBIbJwbQ5vBUfhqbxm0j/4u/n/5alSrLnn76aVJSUvjoo49yfX3cuHEA7N+/n06dOuU6vQe45oVyOHL/PunWrVvRhPWQ7lkUlggxQ8ezMXoId5/7hpVj+pGWnmF1KqWUhQozRfkdd9zBk08+SVxcHO+//z5z5syhVatWNG3alHbt2nH8+HEAJk6cyBNPPAHAgAEDGD16NLfffju1a9e+NF1IUdI/gYuCCDGD3mHTF0HcuesDlo/pRZORUwkJDrI6mVJlz3fPQOKmvNslbnT99GTc4obG8MBrhcuVD9nZ2ZeuzD59+jSdO3dGRBg/fjxvvfUWr7/++lV9jh8/zvLly9m0aRO9e/cu8j0PLRZFRYTGD73K5mlBtN76Five7U6DUdMJCwmxOplSysf06dPn0vODBw/Su3dvEhMTSU9Pp169ern26dq1KyJCTEwMhw8fLvJMWiyKWKPez7FtVgC3bvgPq8d05qYRs4goH2Z1LKXKDk/3AC7uUQz51ntZCigkxx+Zw4cP59lnn6VDhw4sXryY117LffsCAgIuPffG9V86ZuEF9bv9ne1xL9M8Yw37x3bi+KkkqyMppXzU2bNnqVatGsYYPvvsM8tyaLHwkls6jWJP6zeIydpI4vudSHQPSimlSr/CTFF+pRdeeIFu3brRokULKleuXIQp80dKy3QVcXFxprim6s2PPT9+zo1Ln2CnrQ5hQ+dQo1pVqyMpVeps27aN+vXr569TCT4M5S25/Z5EZI0xJi6vvrpn4WV12gzkcLvx1DV7SZnYgX0HD1odSSkFriJRhgpFYWmxKAbRrXtzrMOn1DSHcX7agV17d1sdSSml8kWLRTGp0bILp7pMpgon8Pu8E9t2bLM6klJKeUyLRTGq2rQ953t+RRRnCP2iMxs3bbQ6klJKeUSLRTGr3OgeUvvNJExSiJrelbVrV1sdSSml8qTFwgKRN99O1sNzCbJlUn1OT1auXG51JKXKnCELhjBkwRCrY/gMLRYWqVCnOTJ4Pg6bUHd+H375eYnVkZRShXBxivL169dz22230bBhQ2JiYvjqq6+ualsUU5QDrF27lgULFhRJ/rzodB8WKl+zMcmPLCB1YicaLx7A0sxPuadte6tjKaUKITg4mM8//5y6dety5MgRmjdvzv333094ePilNp5OUZ6XtWvXsnnzZtq39/73hu5ZWKxc1VsIeXQR6Y5Qmv80mO8XzLE6klKqEOrVq0fdunUBqFq1KpUqVeLEiRMe99+1axf3338/zZs356677mLnzp0ATJ06lUaNGhEbG0ubNm1ITU3lpZdeYsqUKQXaK8kvr+5ZiEh74F3ADkw0xlw1A5aI9AZeAAywwRjT370+G7g4z/BBY0xnb2a1UnDl2tge/57THz7A7b89wvyMdDp07m11LKV80uurXmf7qe15trvYxpNxi1sq3MLfW/4931lWrVpFRkYGderU8bjPsGHDmDhxInXq1GH58uWMGDGCRYsW8eKLL7J06VIqV67MmTNnCAoK4rnnnmPz5s2MGTMm39nyy2vFQkTswDjgPiABWC0ic40xW3O0qQv8A2htjDktIpVyvEWqMaaJt/KVNIEVb6TCiO85Oe4B2q55nG8y03iwx0CrYymlCujo0aM8/PDDfPbZZ9hsnh3EOXPmDCtWrKBHjx6X1mVlZQHQunVrBg4cSK9evejevbtXMl+PN/csWgK7jTF7AURkKtAF2JqjzSPAOGPMaQBjTJmebS8gvCqVRy0m8f323L/xCWZnptGlzyOIiNXRlPIZnu4BXNyjmNR+UpFnOHfuHB07duSVV17h1ltv9bifMYbIyMhcxzAmTJjAypUrmTdvHs2aNWPdunVFGTlP3hyzqAYcyrGc4F6XUz2gnogsF5EV7sNWFwWKSLx7fdfcPkBEhrnbxOfnmGBJ5giNosqoxSQG30zHbc8wc/J7XpmbXinlHRkZGXTr1o2BAwfSs2fPfPWNiIigSpUql26L6nQ62bBhAwB79+7l1ltv5eWXXyYiIoLDhw8TGhrK+fPni3wbcmP1ALcDqAvcA/QDJojIxVMGarpnQuwPjBGRqw76GWM+NsbEGWPioqKiiiuz19lDIqg+agFHyjWi657nmP7pmzidWjCU8gXTpk1j2bJl/Pe//710Smx+znaaOnUq48ePJzY2loYNGzJv3jwAnnzySRo3bkzjxo1p06YNjRo1om3btmzYsIGmTZv69AD3YaBGjuXq7nU5JQArjTGZwD4R2YmreKw2xhwGMMbsFZGlQFNgjxfzlii2oPLcOGo+B8Z1pdehV/h6Qjrdhv4fDrvV9V0plZvk5GQABgwYwIABAzzqEx0dzebNmy9bV7t27VzvfzF37tyr1kVFRVFct2bw5jfPaqCuiNQSEX+gL3Dl1s7GtVeBiETiOiy1V0QiRCQgx/rWXD7WUSZIQDlqjpjL/ojW9Dr6/5g1/l9kZjutjqVUqTCp/SSvjFeUVl4rFsaYLGAEsBDYBkwzxmwRkZdE5OJpsAuBJBHZCvwIPG2MSQLqA/EissG9/rWcZ1GVJeIfTPTw2eyPakuvE+8z5/2nScvMtjqWUqqM0Tvl+YrsTPZNGECtxAXMDBvAA8PfJShAL8BXCgp4p7wySO+UVxbY/ag17Av21+hK93P/Y9HYxzifmmF1KqVUGaHFwpfY7EQPmcT+Wn3pcuFrlo79E2cupFudSilVBmix8DU2G9EDx3Og3mAeTJ3Lb+8+zMnzqVanUsrnHHh4IAce1lkSPKXFwheJULPfGA42/DMPZCxk7bt9OXYm2epUSpVpxT1F+axZs3jzzTeLLH9edITUV4lwY6/XOOQXRLv1b/Pje73YhiHKnkzDZ3+xOp1SZVZRTlGelZWFw5H713S3bt2KPvx16J6Fj6vR9XkOt3iWNtm/QkYKR7LKWx1JqTKtsFOU33HHHTz55JPExcXx/vvvM2fOHFq1akXTpk1p164dx4+7ptCbOHEiTzzxBOC6EHD06NHcfvvt1K5d+9J0IUVJ9yxKgWod/85R/yDuWf4vlmc15ExyCuHlgq2OpZQlEv/zH9K35T1Fedp2VxtPxi0C6t/CDc8+m+8sBZmiHCA7O/vSldmnT5+mc+fOiAjjx4/nrbfe4vXXX7+qz/Hjx1m+fDmbNm2id+/eRb7noXsWpUSV+0axXm6htW0L8d9ffYxUKVW8Lk5RPmnSJI+nKL+oT58+l54fPHiQdu3a0bhxY95++222bNmSa5+uXbsiIsTExHD48JUzKxWe7lmUIn4OB6czyhGw9WtM18E6tbkqkzzdA7i4R1Fz8udFnqGgU5RfFBIScun58OHDefbZZ+nQoQOLFy/mtdeuuoccAAEBAZeee+Nia92zKEUaVg0nMSCalhkr2bbvoNVxlCqTCjNFeW7Onj1LtWrVMMbw2WefFUHCgtFigesmKJ7cWtEXVA8PIkCy2P3jZKujKFUmFXaK8iu98MILdOvWjRYtWlC5cuUiTJo/OjcU3r1jVrEzhqOvNuFoRgAN/u83Av3sVidSyusKMjeUNw9DlVQ6N1QhpGSmcCT5COlZpWTaDBHSG/SiGTv4ZdVqq9MoVWLVnPx5mSoUhVXmi8W55CRiliYQ8UvuZxj4ohvvHoQT4dxKPRSllCoaZb5YVLhgY9APTuodyGT7qbzPzfYFtogaJJSPo/nZRSScumB1HKWKRWk5pO4thf39lPli4V+9GufKO2i7wfDpoletjlNkQloOoKYc59el862OopTXBQYGkpSUpAXjGowxJCUlERgYWOD30OssgBZzFrPz3j8QPWM1q+9ZTYsbWlgdqdAqxvUgbfHf8N8yDWfXnthses2FKr2qV69OQkJCvqbVKGsCAwOpXr16gftrsQD8KlcmvF9f7p48hQ8XvErcoBm+f0FbQCgnqt1Hm0NLWLnrCLfdXM3qREp5jZ+fH7Vq1bI6RqlW5g9DXVT50T+Dvx+xc7fzw8EfrI5TJCrfOZjyksL2n6ZZHUUp5eO0WLg5Klak4sMDuX2b4av5b5DlzLI6UqH5123LOUdFbkz4hrOpmVbHUUr5MC0WOUQNHQpBgbT+LoE5u+dYHafwbHbS6/fgLlnPotWbrU6jlPJhWixysIeHU+mPf6LVTsM388eQmuX7tyuNbD0IP8nm9MovrY6ilPJhWiyuUGHQIExYCPd9n8QX276wOk6hyQ2NSCp3My3Pf8+OxPNWx1FK+SgtFlewh4ZSaeijNNtjWPrdR5xNP2t1pEILjOtPE9telvz8s9VRlFI+SotFLioMeAgiyvPgD8l8svkTq+MUWkjzfjixEbB1GhlZTqvjKKV8kBaLXNiCg6n82J9pfMCwfv5kEi8kWh2pcEIrc7rKndzvXMaSbT6+LUopS2ixuIbwvn2RSpF0/ymd8es/tDpOoYXfOoBqksTGX761OopSygdpsbgGW0AAlR8fzs0Jhr2LZrL37F6rIxWKvX4n0u3B1D4yl2Pn0qyOo5TyMVosriO8e3dsVavQd5mT99aMtTpO4fgHk1GvM+1tq5i7erfVaZRSPkaLxXWIvz+Vh4+g1lEnZ3/4no0nNlodqVBCWw6gnKRxbNUMnZ1TKZUvWizyUL5LZxw1b6T/z8KY+Hd8+0u2ZmsuBFXhzpTFrDlw2uo0SikfosUiD+JwUGnkKKodz8axdBXLjyy3OlLB2Wz4Ne3HHbZNfPfbOqvTKKV8iBYLD4R1eAD/ujfRf7mNsavewWl891oF/2b9sYshYNtMLqT7/mSJSqni4dViISLtRWSHiOwWkWeu0aa3iGwVkS0i8kWO9YNEZJf7McibOfMiNhtRo0ZR+WQWlX7Zznf7vrMyTuFE1uVCZCyd+JlvNx21Oo1Sykd4rViIiB0YBzwANAD6iUiDK9rUBf4BtDbGNASecK+vADwPtAJaAs+LSIS3snoi9N57CWjQgH6/2flg9Vgys313yu/gFg/RwHaAVb/9ZHUUpZSP8OaeRUtgtzFmrzEmA5gKdLmizSPAOGPMaQBjzHH3+vuB740xp9yvfQ+092LWPIkIlUaPosKpTOr9lsC0nb57QyFp1JNscVDv2Hz2nki2Oo5Sygd4s1hUAw7lWE5wr8upHlBPRJaLyAoRaZ+PvojIMBGJF5H44rj3bshddxHUpAl9Vzj4dM1HXMi84PXP9IqQimTVvpeu9uVMj99vdRqllA+weoDbAdQF7gH6ARNEJNzTzsaYj40xccaYuKioKC9F/J2IEPXEaMLOZNB8RRKfb/nc65/pLQHN+1NJzpAQ/x1Z2b47YK+UKh7eLBaHgRo5lqu71+WUAMw1xmQaY/YBO3EVD0/6WiLk1lsJbtWKPqscfLFuEkmpSVZHKph67cn0C6NtxhKW7fL+XplSyrd5s1isBuqKSC0R8Qf6AnOvaDMb114FIhKJ67DUXmAh0E5EItwD2+3c60qEqNGjCD6Xwd2rUpiwaYLVcQrGEYCtcQ/a2+OZs3Kn1WmUUiWc14qFMSYLGIHrS34bMM0Ys0VEXhKRzu5mC4EkEdkK/Ag8bYxJMsacAl7GVXBWAy+515UIwc2aEXLnnfRcZWfuxqkknE+wOlKB2Jv0I5AMAnbNIyk53eo4SqkSTHx6+ooc4uLiTHx8fLF9Xuqmzezv1Yvpd/uRMqAjr975arF9dpExhox3Yok/HcLWdlMYemdtqxMppYqZiKwxxsTl1S7fexYiYhORsILFKj2CGjei3L1/oMsqWLplHjtO7bA6Uv6J4N/sIW61b+PHVWt9e94rpZRXeVQsROQLEQkTkRBgM7BVRJ72brSSL2rkKPzSsuge7+Ddte9aHadgYnpjwxB7ahEbE3z/fuNKKe/wdM+igTHmHNAV+A6oBTzstVQ+IvDmeoQ90J4HVmezYecy4hOL7zBYkalQi6zqt9LD8TPTVh+0Oo1SqoTytFj4iYgfrmIx1xiTCegxCyByxAjsmdn0iw/knbW+OYW5o2k/6sgRdm/4hdSMbKvjKKVKIE+LxUfAfiAEWCYiNYFz3grlSwJq16Z85860iU/n0N4NLDm4xOpI+degK06bP+2zl7JwS6LVaZRSJZBHxcIYM9YYU80Y08G4HADaeDmbz4gc/jg2JwyML8e7694ly+ljU38HhSO3dKCb4zdmrN5ndRqlVAnk6QD3aPcAt4jIJyKyFmjr5Ww+w79GDcJ79OD2+GTOH9zL3D1XXntY8klsP8I5h//+JRw6lWJ1HKVUCePpYag/uge42wERuAa3X/NaKh8U+dij2MTGI/HhjFs/jrSsNKsj5c9NfyA7qCI97D/z9RrfvMhQKeU9nhYLcf/sAEw2xmzJsU4BflWqEN63L03iTyMJx/hy+5dWR8ofux/2mF7ca1/HwtXbyHb63kC9Usp7PC0Wa0RkEa5isVBEQgGdqvQKkcMewebvz+NrKzJx00TOZfjYOQCxffEnk2YXfuLXPSetTqOUKkE8LRZ/Ap4BWhhjUgB/YIjXUvkoR1QUFQY8RP34E4QeOcunmz61OlL+VGmCM7Ievf1+YVq8HopSSv3O07OhnLimCf+niPw/4HZjzEavJvNRFf70J2zBwYxeU4kp26ZwPOV43p1KChFssf1oyg62blnPmZQMqxMppUoIT8+Geg0YDWx1P0aJyH+8GcxXOSIiqDBoINFrj1LtaAYfbvjQ6kj5E9Mbg9DJ/MzcDUesTqOUKiE8PQzVAbjPGPOpMeZTXPfD7uS9WL6twuDB2MLCGLWmErN2zWLfWR+6dqF8daTWnfQJ+FWn/1BKXZKfWWdz3u60fFEHKU3sYWFU/OMfqbI+gVsS7by37j2rI+VPbD+qOo/if3QNW47o5IJKKc+LxavAOhH5r4h8BqwBXvFeLN9X4eEB2CMiGLm6It8f+J5NJzZZHclz9R/EOILo5fiFr3WgWymF5wPcXwK3AjOBGcBtxpivvBnM19lCQqg4bBgVNh2i5dEQxqwd4zuTDAaEIvUfpIvfCuav20d6lk4uqFRZd91iISLNLj6AKkCC+1HVvU5dR0S/vjiionh0VRirjq7k1yO/Wh3Jc7F9CXYm0yx9NYu3+tAZXUopr3Dk8fpb13nNoPNDXZctMJCKjz1K1sv/pk2rSoxZO4bbqt6GTbx26/OiU/seTLkb6J+ynE/iO9IxporViZRSFrrut5Yxps11HlooPBDeqxeOqlUYsjyA7UnbWLBvgdWRPGOzIzG9aG3WsnnXHo6eTbU6kVLKQp5eZ9E9l8cfRKSStwP6Opu/P1GPP07gzkM8mFiV99a9R2Z2ptWxPBPbD7vJpqPtN2auPWx1GqWUhfIz3cdE4CH3YwLwd2C5iJT526vmpXyXLvjVvJF+PxsOnz/E9F3TrY7kmcoNoXJjBgb/xrT4Q74zQK+UKnKeFgsHUN8Y08MY0wNogGvMohWuoqGuQ/z8iBo+HMeeBPol1mL8hvGkZPrIPSNi+3JT5k4cp3axat8pq9MopSziabGoYYw5lmP5uHvdKcBHjqlYK6xjR/zr1KHLj6mcTkni862fWx3JM417YcRGH/9fdXJBpcowT4vFUhGZJyKDRGQQMNe9LgQ44714pYfY7USNHIkcOMywxPr8d8t/OZXmA3+ph1ZG6rSll/+vfLfpMOfT9G8DpcoiT4vFcGAS0MT9+AwYboy5YIzRe3F7KLTdfQTUr88ffkgiIz2FCRsnWB3JM7H9iMg8Rkz2Fr7deNTqNEopC3h6BbcBfgGWAD8Ay4yOduab2GxEjRoJCUcZlRjDVzu+4nCyD5xldHMHjH8og0NcA91KqbLH01NnewOrgJ5Ab2CliPT0ZrDSqtw99xAYG8OtCw7hlwXj1o2zOlLe/IORBl1o6/yNrQePsfv4easTKaWKmaeHof4P113yBhljBgItgX95L1bpJSJEjRqFOXacpxObMm/vPHac2mF1rLzF9sU/O4X2jjU6uaBSZZCnxcJmjMk5QVBSPvqqK4TcfjvBcXE0/nYHFQhh7LqxVkfKW83WUL4GQ0NXMmPtYTKz9RbsSpUlnn7hLxCRhSIyWEQGA98C870Xq3QTEaJGj8J5Mom/HY5lWcIy1hxbY3Ws67PZIKY3DdPWIMmJLN1xwupESqli5OkA99PAx0CM+/GxMUYvxiuE4BYtCGndmjrfbKCGLZJ31rxT8q+QjumLGCf9g1fpQLdSZYzHh5KMMTOMMX9xP2Z5M1RZETV6FM7TZ/jbgYZsOLGBHw/9aHWk64uqB1Wb0T/wV5ZsP87x82lWJ1JKFZO87mdxXkTO5fI4LyLn8npzEWkvIjtEZLeIPJPL64NF5ISIrHc/huZ4LTvH+rkF27ySLSgmhnJt2lBl7ipu8avB2LVjyXaW8BsNxfajcsou6poDzF7nA6f9KqWKRF5TlIcaY8JyeYQaY8Ku11dE7MA44AFcc0n1E5EGuTT9yhjTxP2YmGN9ao71nfO7Yb4iatRInOfO89fd9dhzdg9z95TwutioB9gc/DliNdPiE0r+oTOlVJHw5hlNLYHdxpi9xpgMYCrQxYuf55MC69cntH17Imb/TIvAW/hgwwekZ6dbHevaQipC3Xa0y/6JfcfPsu6QzvaiVFngzWJRDcg5CprgXnelHiKyUUSmi0iNHOsDRSReRFaISNfcPkBEhrnbxJ844btn50SNHIEzLY1RW2uQeCGRqdunWh3p+mL7EpR+kjZ+W/laB7qVKhOsvlbiGyDaGBMDfI9rzqmLahpj4oD+wBgRqXNlZ2PMx8aYOGNMXFRUVPEk9oKAOnWFwQ0tAAAdv0lEQVQo/2AnAuf8yL3l4piwaQLnM0rwVdL12kNgef5cIZ5vNhwlJSPL6kRKKS/zZrE4DOTcU6juXneJMSbJGHPxmMtEoHmO1w67f+4FlgJNvZjVcpGPP47JzGTY+oqcTT/LpM2TrI50bY4AaNidJsk/Y9LP892mRKsTKaW8zJvFYjVQV0RqiYg/0BfX1OaXiEiVHIudgW3u9REiEuB+Hgm0BrZ6Mavl/GvWJLx7N2T2InqG3cPkrZM5nnI8745Wie2HPTuNh8PW6zUXSpUBXisWxpgsYASwEFcRmGaM2SIiL4nIxbObRonIFhHZAIwCBrvX1wfi3et/BF4zxpTqYgEQ+ec/A9B/lT9ZzizGbxhvcaLrqNESImoxIHgFK/edYv/JC1YnUkp5kVfHLIwx840x9YwxdYwxr7jXPWeMmet+/g9jTENjTKwxpo0xZrt7/a/GmMbu9Y2NMZ94M2dJ4Ve1KuG9e5P1zSIGht/PzF0z2X92v9WxcicCsX2pdiaeapLE9DU6uaBSpZnVA9zqChUfHYbY7XT5OR1/uz/vrXvP6kjXFtMHwfBE5fVMX5NAtlOvuVCqtNJiUcL4VapERP/+pH+7iMfKP8iiA4vYcnKL1bFyV6EW3Hgb7bOXkngulZ93+e7py0qp69NiUQJVfGQoEhhI28UniAiI4J2171gd6dpi+hB6fg+tgw7pfS6UKsW0WJRAjgoVqDDwYVIXfM+osC6sPLqSX4/8anWs3DXsCvYARkWtYdHWRE5dyLA6kVLKC7RYlFAVhwzBFhpKi3l7qBpSlTFrxuA0JfCGQ0ERcHN7mp/7AZOdyZz1OrmgUqWRFosSyl6+PBX/OISUJT/yl3Jd2HZqG4v2L7I6Vu5i++FIO8XAqN18tfqQTi6oVCmkxaIEi3h4IPbwcG6Zvo6bwm/ivXXvkenMtDrW1W66F4IrMjB4BdsTz7PlSJ6z1yulfIwWixLMXi6Eio8MJeWXX3javxMHzx9k5s6ZVse6mt0PGvWk5smfiHSk6BXdSpVCWixKuIj+/bFHRlL1i59oFtWUDzd8SEpmitWxrhbbF8lO56lq25i97jBpmSX8Jk5KqXzRYlHC2YKCiHz0UVJXr+Yv0o6ktCQmb51sdayrVW0KkTfzgPMnzqVlsWjrMasTKaWKkBYLHxDepzeOKlUo/9m3tKl+D5O2TOJ02mmrY11OBGL7UP5EPC3Dzuh9LpQqZbRY+ACbvz+Rjz1G6vr1jExvTWpWKhM2TbA61tUa9waEJyuv55fdJ0k4XQIPlymlCkSLhY8I794Nvxo1sE+cRufaDzJ1+1SOJB+xOtblwmtArTuJO7cIYwwz1ug1F0qVFlosfIT4+RE5/HHSt27jkdMxCMK49eOsjnW1mL74nd3P4BrH+HrNIZw6uaBSpYIWCx9S/sEH8a9dm6yPJ9OvXh++2fMNu07vsjrW5Rp0BkcQg0JWkHA6lRV7k6xOpJQqAlosfIjY7USNHEH6rt08dCSaEL8Qxq4da3WsywWEQv0HiU5cSIVAo9dcKFVKaLHwMaH330/AzTeTMv5T/lh/MEsTlrL22FoAhiwYwpAFQyxOCMT2QdLO8lT0Pr7bnMjZ1BJ41blSKl+0WPgYsdmIGjWSjAMH6LKrPJFBkYxZO6ZkzcdU6x4odwMdzTLSs5x8s6GEDcQrpfJNi4UPKte2LYGNGnFu/AQeb/AI646v46eEn6yO9Tu7A2J6EXZoCS0rOfWaC6VKAS0WPkhEiBo9mswjR2iz0VAzrCbvrn23ZO1dxPRFnFk8UWUzGxLOsiPxvNWJlFKFoMXCR4Xc0Zqg5s05Pf5jRjV4lN1ndpOUVoLOPLqhEVRuTIuzi/Czi+5dKOXjtFj4KNfexSiyjh+n+W+naFixIUeSj5SsGyTF9sUvcS0P1Uln1rrDZGSVoGxKqXzRYuHDQlq2JPi2W0maMIEnGjxGhjOD4ynHrY71u8Y9QWwMLreSpAsZLNlegrIppfJFi4WPixo1iuykJG5avIsw/zCOXjjKlpNbrI7lEnoD1GlLzcPfcEOonx6KUsqHabHwccFNm1Lu7rtJ+uQTbnJUxSY2+s/vzxur3ygZ972I6YucTWD0TSf4ccdxjp1LszqRUqoAtFiUApGjRuI8e5Y7lp+hYcWG9Kzbk8lbJ9NtTjd+TvjZ2nC3dAT/cnTiJ5wGZq7VyQWV8kVaLEqBoIYNCW3XjtsXJTBw7Hb+ddu/+Kz9ZwQ6Ann8h8f5209/42TqSWvC+QdDg66E7pnPHTWD+Tr+UMk6xVcp5REtFqVE1MgRiBNCz2QA0KxyM75+8Gsej32cxQcX02V2F2btmmXNF3VsH8g4z4iqO9l78gJrDpSwGzcppfKkxaKUCKhbl5RyDkLPZnBhxUoA/O3+/LnJn5n+4HRuCr+J5359jqGLhnLg3IHiDVfzDgirTtzZhYT423VyQaV8kBaLUuRMhQCyHTYODhnCibHvYbKyAKgdXptJ7Sfx3G3PsS1pG93ndGfCxglkOotpgj+bDWL74Nj3I33rBzBv41EupGcVz2crpYqEFotSJNvPRmL1YMp37crJDz7gwODBZB49CoBNbPSq14s5Xedwd427GbtuLL2/6c3GExuLJ1xMXzBOBoetJiUjm283HS2ez1VKFQktFqXI1JEN+XJ0I6q++h+qvvE66Vu3sa9rN84vWXKpTVRwFG/f8zZj24zlXMY5BswfwKsrX+VC5gXvhouqB1WbUf3gHGpHheg1F0r5GC0Wpcik9pOY1H4SAOU7d6bWzBn4VatGwuPDSXzlPzgzMi61bXNjG+Z0mUPfW/ry5fYv6TK7C0sPLfVuwNh+yLHNPHZLKqv3n2bviWTvfp5Sqsh4tViISHsR2SEiu0XkmVxeHywiJ0RkvfsxNMdrg0Rkl/sxyJs5Syv/6GhqTv2SCoMGcnryZPb37Uv6vn2XXi/nX45nWz3L5A6TCfUPZeSSkfx16V85kXLCO4Ea9QCbg47OZdhtwtdrErzzOUqpIue1YiEidmAc8ADQAOgnIg1yafqVMaaJ+zHR3bcC8DzQCmgJPC8iEd7KWprZ/P2p/I9/UP2DD8g6cpR9PXpyds6cy9rERsUyrdM0RjUdxdJDS+kyuwvTd04v+kkJQypC3XaE7JhJ23oVmLEmgaxsnVxQKV/gzT2LlsBuY8xeY0wGMBXo4mHf+4HvjTGnjDGnge+B9l7KWSaEtm1DrdmzCGrQgCN/f4Yjf38G54Xfxyn87H48EvMIMzrP4JaKt/Diby8yZMEQ9p7dW7RBYvtCciKP1kjg+Pl0lu3y0l6MUqpIebNYVANyjmImuNddqYeIbBSR6SJSIz99RWSYiMSLSPyJE/qlkxe/G27gxs/+S+Tw4Zz95hv29ehJ2tatl7WJLh/NJ+0+4aXbX2L3md30nNuTDzd8SGZ2EZ1mW689BJan2ekFRJbzZ9pqPRSllC+weoD7GyDaGBODa+/hs/x0NsZ8bIyJM8bERUVFeSVgaSN2O1EjR3DjpEk4U1LY36cvpyb/77Iru0WEbnW7MafrHO698V4+WP8Bvb7pxbrj6wofwBEADbtj2/EtfWLCWbztGEnJ6YV/X6WUV3mzWBwGauRYru5ed4kxJskYc/GbYiLQ3NO+qnBCWrWk1pzZhNx+O8deeYWEESPJOn35NByRQZG8cfcbjPvDOFKyUhj43UD+veLfnM8o5C1SY/tBZgoPl99EltMwa53+0ypV0nmzWKwG6opILRHxB/oCc3M2EJEqORY7A9vczxcC7UQkwj2w3c69ThUhR0QE1cd/SOV/PEPysmXs69adlPj4q9rdVf0uZneZzYD6A/h659d0nd2VHw78UPAPrtESImpxw/5ZNKkRzjSdXFCpEs9rxcIYkwWMwPUlvw2YZozZIiIviUhnd7NRIrJFRDYAo4DB7r6ngJdxFZzVwEvudaqIiQgVBg0i+ssvEX9/DgwcxIkPPsBkZ1/WLtgvmL+3/DtTOkwhIjCCJ5Y+weglozl24VhBPtQ10L3vZwY39GPnsWQ2Jpwtoi1SSnmDlJa/6OLi4kx8Ln8VK89lJyeT+MKLnJs3j+BWraj6xhv4Va50VbtMZyaTt07mg/Uf4LA5eKLZE/S+uTc2ycffHqf2wtimpN39L5osaUiPZtV5pVvjItwapZQnRGSNMSYur3ZWD3CrEsRerhxV33yDKv/5D6kbN7Kva1eSf/rpqnZ+Nj/+2OiPzOo8i8aRjXll5SsM+m4Qu0/v9vzDKtSGGrcSuGUaDzS8gbnrj5CakZ13P2WpVpN60GpSD6tjFFpp2Q4ovm3RYqEuIyKEd+9GrRnTcVSqxKFHH+PY629gckwVclGNsBp8fN/HvHLHK+w/t59e83rx/rr3Sc/28Oym2L5wcgeDa5/hfHoWC7ckFvHWlAxDFgxhyIIhVsdQqlC0WKhcBdSuTfS0r4jo359Tkyaxv/9DZBw8eFU7EaFznc7M6TqH9tHt+WjjR/Sc25P4RA8OCTbsCnZ/Yk5+R40KQZfd56I0fcFuPXqOrUfPWR1DqULRYqGuyRYQwA3P/Ytq740l4+BB9nXrztlvv821bYXACrx656uMv3c8mc5Mhiwcwgu/vsDZ9OsMXAdFwM0PIJtn0KdpFX7dk8ShUyle2hqlVGFosVB5CrvvPmrPmklAvXoc+etTHPnnP3Gm5P6l3rpaa2Z2nsmQhkOYvXs2XWZ3YeH+hdc+NTa2H6ScpF/FnYhwaXJB/WtcqZJFi4XyiF+1atSc/DkVH3uUszNmsq9Xb9J27My1bbBfMH+J+wtfdvySSsGVeOqnpxi5ZCSJF3IZk7jpXgiuSMXdM7njpkimxx8i21k6ztBTqjTRYqE8Jg4HlZ54ghs//YTsc2fZ37s3p6dOveZeQ/2K9fmi4xc8FfcUqxJX0WV2F6Zsm0K2M8dZT3Y/aNQTdnxH/5jyHDmbxq97ThbTFimlPKXFQuVbyG23UXv2bIJbtCDxhRc5PPoJss/lfsjIYXMwqOEgZnaeSdNKTXlt1Ws8/N3D7Di14/dGsX0hO517za+UD/JjWrxOLqhUSaPFQhWIo2JFanz8EZWeforzS5awr2s3Utevv2b76qHV+fDeD3ntztc4nHyYvvP68u7ad0nLSoOqTSGyHn6bp9G1SVUWbknEme1fjFujlMqLFgtVYGKzUfFPfyJ6yv9AhP0PDeDkhAkYZ+43NBIROtbuyJwuc+hYuyMTN02kx9werExc5dq7OPgb/es5ychyknauTjFvjVLqerRYqEILio2l1uxZhN53HyfeeptDQx8h6+S1xx3CA8P59x3/ZmK7iRgMQxcN5Z+ZCZyx2bn52Hc0qBJG6tm6xbgFSqm8aLFQRcIeGkq1d97mhpdeJGXNGvZ27Uby8uXX7dOqSitmdp7J0MZD+TZhKV1urMG3W6fQq1lVstIjyUyrUEzplVJ50WKhioyIENG7N9FfT8MREc6hPw3l+FtvYzKvfZe9QEcgo5uNZmqnqVQLrswzIYZfT/0f4n+SC6diWLbzBCv2JrHu4Gm2HjnH7uPJHDqVwvFzaZxNySQtMxunnmqrlNc5rA6gSp/AevWInjaNY6++RtKECaSsWkXVt97Cv3pud9V1ubnCzUx+cBpTP2rKu2wnpNYOMs/EMfSbHYAANjDiem5sl5aNe9lus+Gw2bHb7PjZ7Djcy352Bw6bDT+7a72f3bXOz2bH3+5a9nc48Lc78Le71gU4HATYHfj72fG3Owh0OPC3+xHgZyfQ4SDA4SDQz7U+0OEgIMdzu82OTWzYxfXTJjYMrmKW6cwEA+bi/4zBaVzjOxeXc7521bL7fS72+/09DFlOJ07jxGkg25lNtnt9ttOJ0xicTkO2cbXJdj83xrjaOS+2c17q57z4+qV+rnXpaYEATF73Yy7/iiaXZzlezeUUa/P7i9f8b+N672py+ySTx+tAmns7Jsb7/m1y0tICEfH+H0w6RbnyqnPffcfRfz0HIlT5978Ju7/d9TvMGMrRPYvpEFWFLEnBYbPjxIlxfzkqpa7mTLuBLY9+X6C+nk5RrnsWyqvCHniAwMaNOfyXv3J49Ggu9O1D5WeewRYYmHuH2L5U2fQ17VJqsDT4JlYOnHHppYt/Ubv+As7O9fnFZddfzdm5LufWL8uZTXpWFhnZrp9pWZlkZmWTnp1NxsX12Vlkup9nZmeTke16nuXMvrQu0+n6meXMIis7myzjZNepAwDUCq+ODUHEhoi4nwvi/mlzP+Dicxvifn7xdXH3t4sgAoLtstddezNyWT/7pddt7geIu50tZz+bK5Or3++v222/Px8X/z8ARrUYmOs/nyA5F67/+hUNc3vl9ybXbiNy9Vq5zhLAWysmAvD0bY9c71N9wpu/TcBRDHsWWiyU1/lXr070/yZz/N13OfXJp6SuWUu1d94m4Kabrm5c6x4oV5kOyUdZGnz5jZdcX3x27Njxw694wheBi/camD/wU4uTFN4nm11fsn+Ku9/iJIXz/vr3AXioyT3WBikCY9e9VyyfowPcqliIvz+Vn36aGhM+JispiX09e3Fm+vSrj2PbHdC4F61TT9AoLfe5p5RSxU+LhSpW5e68k1qzZxHUtAlH//kvjvz1KbKTky9vFNsPB3B72gVLMiqlrqaHoVSx86tUiRsnTiRpwkROvPceqZs2Ue3ttwhq7L4H9w2N2O/wp/e50zClF4TXhIhoiHD/DK8JgWFWboJSZY4WC2UJsduJfOxRglu24PBTT7G/X38q/eUvVBg8CLHZ+Kx8BdpfOEer80fh4ApIv2KiwqAKruKRWyEpXwMcOreUUkVJi4WyVHCzZtSeOZOj//oXx994gwsrV1D11VdpMcdOEhGw+BfXOfipp+HMATi9H067f545AImbYPu34Mxx4Z/YIKxa7oUkIhrKVbp0do1SyjNaLJTl7OHhVBs7ltNffsnx115nX5euBBpDWqD7C10Egiu4HlWbXv0Gzmw4f/T3QpKzqOxeDMlX3HTJEQThN+ZeSCJqQkCoV7dXKV+kxUKVCCJChf79CW7WjMNP/oUb9p0gpZydY2+8iS042PUICXH/vGI5OBhbcCi2aq2Q6DuufvPMVDhz8IpCst+1fPC3axziis69kJSv4bphUz48P2Wb68mQfP9aSpzSsi2lZTug+LZFi4UqUQJvuYVaM6az5u5bCbqQxekvvsCkpXncX/z9cxSXYGzBIdhCgpHLCk44tuDbsQXfh612EDY/sGWfx5Z9BltmErb0E9jSE7HtX49tyzyEKw9xVXcXkpoQHv17IQmvqYe4VKmlxUKVOLbgYE5Vcl3hff+3qzDZ2ThTU3FeuIDzQgrOlBScKRdcPy+4fpqUlMuWf2/nXncy6ffllJR8FKAoxN8PW1AANn87Nj/BZs/GZkvAxm5skoo4DDaHweZwYgvwwxYWgS0iClvEDdgiq3FDYBpZxkbaD1O4dDXxpYJyjWVbjuVLr9kua3LpymfbFX1zvuf13je3z8/Z3nbFmfUihNld81ulr/3p2r8yHxBmc/0BkLF+mcVJCi/UlkV23s0KTeeGUiXSgYdd00nUnPy5V97fZGW5CtDFwnKxyFwsQh4XomScyck4U1MxGdeeXVcpb/ILz+SmFbsL1FfnhlI+zVtF4iJxOLCHhmIPLbrB7EsF6EKOvZ7Tx9jy4p+xiaHekFHuhhd7OHMsu1de+uPtiuXL/qgzv//IbT0GnFe8z6X3yNE3xyy2V70v7rZX/C15YPYUAGp26Z/Lb8B3HJjzBQA1O/ezOEnhHZj7JSnF8DlaLJQqIrkXoMactrmu+Qh7aJQ1wYpQwozpAMQ8/rLFSQonYaZrgsqY4f+2OEnhJcyaWSyfo8VCKS+7hdJzgWBDCbA6QpEoLdsBxbctOmahlFJlmKdjFjqRoFJKqTxpsVBKKZUnrxYLEWkvIjtEZLeIPHOddj1ExIhInHs5WkRSRWS9+zHemzmVUkpdn9cGuEXEDowD7gMSgNUiMtcYs/WKdqHAaGDlFW+xxxjTxFv5lFJKec6bexYtgd3GmL3GmAxgKtAll3YvA68Dns/poJRSqlh5s1hUAw7lWE5wr7tERJoBNYwx3+bSv5aIrBORn0Tkztw+QESGiUi8iMSfOHGiyIIrpZS6nGUD3CJiA94G/prLy0eBG40xTYG/AF+IyFW3RjPGfGyMiTPGxEVFRXk3sFJKlWHeLBaHgRo5lqu7110UCjQClorIfuBWYK6IxBlj0o0xSQDGmDXAHqCeF7MqpZS6Dm8Wi9VAXRGpJSL+QF9g7sUXjTFnjTGRxphoY0w0sALobIyJF5Eo9wA5IlIbqAvs9WJWpZRS1+G1s6GMMVkiMgJYCNiBT40xW0TkJSDeGDP3Ot3vAl4SkUxcs609Zow5db3PW7NmzUkROVCIyJHAyUL0LylKy3aAbktJVVq2pbRsBxRuW2p60qjUTPdRWCIS78kl7yVdadkO0G0pqUrLtpSW7YDi2Ra9glsppVSetFgopZTKkxaL331sdYAiUlq2A3RbSqrSsi2lZTugGLZFxyyUUkrlSfcslFJK5UmLhZuIvCwiG92z3C4SkapWZyooEXlTRLa7t2eWiIRbnamgRKSXiGwREefFWYl9iaczL/sCEflURI6LyGarsxSGiNQQkR9FZKv7v63RVmcqKBEJFJFVIrLBvS0veu2z9DCUi4iEGWPOuZ+PAhoYYx6zOFaBiEg7YIn7WpfXAYwxf7c4VoGISH1c19p8BDxljPGZ2yG6LyzdSY6Zl4F+V8687CtE5C4gGfjcGNPI6jwFJSJVgCrGmLXuWa/XAF198d9FRAQIMcYki4gf8Asw2hizoqg/S/cs3C4WCrcQwGerqDFmkTEmy724AtdUKz7JGLPNGLPD6hwF5OnMyz7BGLMMuO7Fsb7AGHPUGLPW/fw8sI0rJjn1FcYl2b3o53545btLi0UOIvKKiBwCHgKeszpPEfkj8J3VIcqoPGdeVtYSkWigKVffT8dniIhdRNYDx4HvjTFe2ZYyVSxEZLGIbM7l0QXAGPN/xpgawBRghLVpry+vbXG3+T8gC9f2lFiebItSRU1EygEzgCeuOLLgU4wx2e4bxVUHWoqIVw4Rem1uqJLIGHOvh02nAPOB570Yp1Dy2hYRGQx0Av5gSvjAVD7+XXxNXjMvK4u4j+/PAKYYY2ZanacoGGPOiMiPQHugyE9CKFN7FtcjInVzLHYBtluVpbBEpD3wN1yz+KZYnacMu+7My8oa7kHhT4Btxpi3rc5TGO4ZusPdz4NwnUzhle8uPRvKTURmADfjOvPmAK6Zbn3yr0AR2Q0EAEnuVSt8+MyubsB7QBRwBlhvjLnf2lSeE5EOwBh+n3n5FYsjFZiIfAncg2uG02PA88aYTywNVQAicgfwM7AJ1//fAZ41xsy3LlXBiEgM8Bmu/75swDRjzEte+SwtFkoppfKih6GUUkrlSYuFUkqpPGmxUEoplSctFkoppfKkxUIppVSetFgolQ8ikpx3q+v2ny4itd3Py4nIRyKyR0TWiMhSEWklIv4iskxEytRFs6pk02KhVDERkYaA3Riz171qIq6J+eoaY5oDQ4BI96SDPwB9rEmq1NW0WChVAOLypnsOq00i0se93iYiH7jvJ/K9iMwXkZ7ubg8Bc9zt6gCtgH8aY5wAxph9xphv3W1nu9srVSLobq5SBdMdaALE4rqiebWILANaA9FAA6ASrumvP3X3aQ186X7eENfV6NnXeP/NQAuvJFeqAHTPQqmCuQP40j3j5zHgJ1xf7ncAXxtjnMaYRODHHH2qACc8eXN3Eclw35xHKctpsVCq+KQCge7nW4BY9930riUASPN6KqU8oMVCqYL5GejjvvFMFHAXsApYDvRwj11UxjXx3kXbgJsAjDF7gHjgRfcsqIhItIh0dD+vCJw0xmQW1wYpdT1aLJQqmFnARmADsAT4m/uw0wxcd8TbCvwPWAucdff5lsuLx1CgMrBbRDYD/8V1tzOANu72SpUIOuusUkVMRMoZY5LdewergNbGmET3/QZ+dC9fa2D74nvMBJ4xxuwshshK5UnPhlKq6M1z35DGH3jZvceBMSZVRJ7HdR/ug9fq7L5R0mwtFKok0T0LpZRSedIxC6WUUnnSYqGUUipPWiyUUkrlSYuFUkqpPGmxUEoplSctFkoppfL0/wFSB7xyOGGk7wAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# plot CV误差曲线\n",
    "test_means = grid.cv_results_[ 'mean_test_score' ]\n",
    "test_stds = grid.cv_results_[ 'std_test_score' ]\n",
    "train_means = grid.cv_results_[ 'mean_train_score' ]\n",
    "train_stds = grid.cv_results_[ 'std_train_score' ]\n",
    "\n",
    "\n",
    "# plot results\n",
    "n_Cs = len(Cs)\n",
    "number_penaltys = len(penaltys)\n",
    "test_scores = np.array(test_means).reshape(n_Cs,number_penaltys)\n",
    "train_scores = np.array(train_means).reshape(n_Cs,number_penaltys)\n",
    "test_stds = np.array(test_stds).reshape(n_Cs,number_penaltys)\n",
    "train_stds = np.array(train_stds).reshape(n_Cs,number_penaltys)\n",
    "\n",
    "x_axis = np.log10(Cs)\n",
    "for i, value in enumerate(penaltys):\n",
    "    #plt.plot(log(Cs), test_scores[i], label= 'penalty:'   + str(value))\n",
    "    plt.errorbar(x_axis, -test_scores[:,i], yerr=test_stds[:,i] ,label = penaltys[i] +' Test')\n",
    "    plt.errorbar(x_axis, -train_scores[:,i], yerr=train_stds[:,i] ,label = penaltys[i] +' Train')\n",
    "    \n",
    "plt.legend()\n",
    "plt.xlabel( 'log(C)' )                                                                                                      \n",
    "plt.ylabel( 'logloss' )\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "上图给出了L1正则和L2正则下、不同正则参数C对应的模型在训练集上测试集上的logloss。可以看出在训练集上C越大（正则越少）的模型性能越好；但在测试集上当C=1 时性能最好（L1正则）"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [],
   "source": [
    "lr_best = grid.best_estimator_"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 用LogisticRegressionCV实现正则化的 Logistic Regression\n",
    "\n",
    "\n",
    "## L1正则 Solver='libliner' scoring = 'neg_log_loss'"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "LogisticRegressionCV(Cs=[0.001, 0.01, 0.1, 1, 10, 100, 1000],\n",
       "           class_weight=None, cv=5, dual=False, fit_intercept=True,\n",
       "           intercept_scaling=1.0, max_iter=100, multi_class='ovr',\n",
       "           n_jobs=None, penalty='l1', random_state=None, refit=True,\n",
       "           scoring='neg_log_loss', solver='liblinear', tol=0.0001,\n",
       "           verbose=0)"
      ]
     },
     "execution_count": 13,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "from sklearn.linear_model import LogisticRegressionCV\n",
    "\n",
    "Cs = [1e-3, 1e-2, 1e-1, 1, 10, 100, 1000]\n",
    "#nCs = 1  #Cs values are chosen in a logarithmic scale between 1e-4 and 1e4.\n",
    "\n",
    "# LogisticRegressionCV比GridSearchCV快\n",
    "lrcv_L1 = LogisticRegressionCV(Cs=Cs, cv = 5, scoring='neg_log_loss', penalty='l1', solver='liblinear', multi_class='ovr')\n",
    "lrcv_L1.fit(X_train, y_train)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "{1: array([[-0.69314718, -0.6359066 , -0.48852103, -0.4876455 , -0.48865214,\n",
       "         -0.48877697, -0.48878928],\n",
       "        [-0.69314718, -0.64460736, -0.52488541, -0.52641178, -0.52786951,\n",
       "         -0.52802904, -0.52804445],\n",
       "        [-0.69314718, -0.63328109, -0.45813512, -0.45606985, -0.45653936,\n",
       "         -0.45659911, -0.45660557],\n",
       "        [-0.69314718, -0.63024185, -0.43303243, -0.42140792, -0.42098201,\n",
       "         -0.42094873, -0.4209422 ],\n",
       "        [-0.69314718, -0.63580105, -0.4865079 , -0.48335898, -0.4827114 ,\n",
       "         -0.48265649, -0.48265967]])}"
      ]
     },
     "execution_count": 14,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "lrcv_L1.scores_"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "验证只有Target=1的值， 0的值不考虑？"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYsAAAEKCAYAAADjDHn2AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzt3XucVPV9//HXZ2Z3WWBhuewCu4sFVEBQYIgr+gu5mNYINcpuG5tIrNUYtX0k1uTXR9OYJr+fv2japLc0aX/20SjamIua/Ey7YDQaGjEmaUEWXRQWQcSI3Jflftn75/fHHHRE2Jm9nD1zeT8fj3kw5zvnzLyP4rw9lznH3B0REZHexKIOICIi2U9lISIiaaksREQkLZWFiIikpbIQEZG0VBYiIpKWykJERNJSWYiISFoqCxERSasozDc3s8XAt4A4sMzdv37a6/8IfCiYHAFMcPcxwWs3Al8OXvuquz/U22dVVFT41KlTBzG9iEj+W7du3X53r0w3n4V1uQ8ziwNbgA8DO4C1wFJ3bz7L/H8KzHf3m81sHNAI1AIOrAMudveDZ/u82tpab2xsHOS1EBHJb2a2zt1r080X5m6oBcBWd9/m7h3Ao0BdL/MvBR4Jni8CVrr7gaAgVgKLQ8wqIiK9CLMsaoA3U6Z3BGPvYmZTgGnAM31dVkREwpctB7ivAx5z9+6+LGRmt5lZo5k1trS0hBRNRETCLIudwDkp05ODsTO5jrd3QWW8rLvf5+617l5bWZn2+IyIiPRTmGWxFphuZtPMrIRkIaw4fSYzuwAYC/x3yvDTwJVmNtbMxgJXBmMiIhKB0E6ddfcuM7ud5Jd8HHjQ3Tea2d1Ao7ufKo7rgEc95bQsdz9gZveQLByAu939QFhZRUSkd6GdOjvUdOqsiEjfZcOpszmhu8f52pOb2HnoZNRRRESyVsGXxfYDJ3j4+e0svW81uw+rMEREzqTgy2JaxUi+e/MCDhzv4BP3r2HfkbaoI4mIZJ2CLwuA+b81loduvoR9R9pYev9qWo62Rx1JRCSrqCwCF08Zx4M3XcKuQ21cv2w1rcdUGCIip6gsUlx67ngeuLGWN1pPcP2yNRw83hF1JBGRrKCyOM17z69g2Y21bNt/nBseXMPhE51RRxIRiZzK4gzeP72Sb99wMVv2HOOPHlzDkTYVhogUNpXFWXxo5gT+5fr3sHHXEW568HmOtXdFHUlEJDIqi15cMXsi//cT81m/4zA3/9taTnSoMESkMKks0lh8URXfui5B4xsHuPk7aznZ0aerqIuI5AWVRQaunlvNNz6WYM3rB7j1u420daowRKSwqCwyVD+/hr+7dh6/fm0/f/y9dbR3qTBEpHCoLPrg2osn87Xfm8MvtrTw6e+/QEdXT9SRRESGhMqij65b8FvcU38RP39lH3/6yAt0dqswRCT/qSz64YbLpvB/rpnN0xv38rlHm+hSYYhIngvtTnn57qaF0+jqcb76xCaK4sY3PpYgHrOoY4mIhEJlMQC3vP9cOrp7+NunNhOPGX9/7TxiKgwRyUMqiwH69OXn09XtfGPlFopjMb72+3NUGCKSd1QWg+CO35lOV3cP//TMVorixlfrL8JMhSEi+UNlMUj+54dn0NHt/OsvXqM4HuOua2arMEQkb6gsBomZ8YXFM+nq7mHZr16nKGZ86SOzVBgikhdUFoPILFkQXT2eLIx4jC8snqnCEJGcp7IYZGbGXdfMprO7h3/9xWuUxI0/u3Jm1LFERAZEZRECM+Oeuovo7vHgoHeMO35netSxRET6TWURkljM+Ovfm0NncFptUdz49OXnRx1LRKRfVBYhisWMv712Ll09yR/ulcRj3PL+c6OOJSLSZyqLkMVjxj/8wby3Lg0SjxmfXDgt6lgiIn2ishgCRfEY3/x4gq7uHr7yeDNF8Rg3XDYl6lgiIhnTVWeHSHE8xj8vfQ9XzJrA/2rYwKPPb486kohIxlQWQ6ikKMa917+HD86o5Iv/8TKPrdsRdSQRkYyoLIbYsKI4377hYhaeV8HnH1tPw4s7o44kIpKWyiICpcVx7v+jWi6bNp4/+1ETP3lpV9SRRER6pbKIyPCSOA/cVEvtlHF89tEmntqwO+pIIiJnpbKI0IiSIh785CXMm1zO7Q+/yH827406kojIGaksIlY2rIjv3LyAC2vK+fQPXmDV5n1RRxIReReVRRYYXVrMdz+5gBmTyvjj763juS0tUUcSEXmHUMvCzBab2WYz22pmd55lno+ZWbOZbTSzh1PGu82sKXisCDNnNigfUcz3br6U8yrLuPW7jfzX1v1RRxIReUtoZWFmceBe4HeB2cBSM5t92jzTgS8CC939QuBzKS+fdPdE8FgSVs5sMnZkCd//1AKmjB/Bpx5qZM221qgjiYgA4W5ZLAC2uvs2d+8AHgXqTpvnVuBedz8I4O4Fv8N+fNkwfnDLZVSPKeWT31nLujcORB1JRCTUsqgB3kyZ3hGMpZoBzDCzX5vZajNbnPJaqZk1BuP1Z/oAM7stmKexpSV/9vNXjhrGI7dexsTRpdz44Fqa3jwUdSQRKXBRH+AuAqYDlwNLgfvNbEzw2hR3rwU+AXzTzM47fWF3v8/da929trKycqgyD4kJo0t5+NZLGTeyhBseWMPLOw5HHUlECliYZbETOCdlenIwlmoHsMLdO939dWALyfLA3XcGf24DngXmh5g1K1WVD+eR2y6jfHgxf/jAGpp3HYk6kogUqDDLYi0w3cymmVkJcB1w+llNDSS3KjCzCpK7pbaZ2VgzG5YyvhBoDjFr1qoZM5xHbr2MkSVx/vCBNWzeczTqSCJSgEIrC3fvAm4HngY2AT9y941mdreZnTq76Wmg1cyagVXA5929FZgFNJrZ+mD86+5ekGUBcM64ETx862UUx43rl61m6z4VhogMLXP3qDMMitraWm9sbIw6RqheaznGdfetBuCHt13GuZVlEScSkVxnZuuC48O9ivoAt/TBeZVlPHzLpfT0OJ+4fw1vtB6POpKIFAiVRY6ZPnEUD996Ge1d3Sy9bzVvHjgRdSQRKQAqixw0c9Iovn/LpRzv6Gbp/avZeehk1JFEJM+pLHLUhdXlfP9Tl3L4ZCefuH81ew63RR1JRPKYyiKHzZlczndvXkDrsQ6W3r+afUdUGCISDpVFjpv/W2N56OZL2HukjaX3r6blaHvUkUQkD6ks8sDFU8bxbzddwq5DbfzhsjWc7OiOOpKI5BmVRZ649NzxfOu6BJv3HuVnzXuijiMieUZlkUeumDWR6vJSljftijqKiOQZlUUeicWMaxLVPLelhQPHO6KOIyJ5RGWRZ+oTNXT1OE+8pK0LERk8Kos8c8GkUcyYWEaDdkWJyCBSWeQZM6MuUcO6Nw7qUiAiMmhUFnmoLlENwPKm0+81JSLSPyqLPDR57AgumTqWhqZd5Msl6EUkWiqLPFWXqGHrvmM079atWEVk4FQWeeojc6ooipl+cyEig0JlkafGjizh8pmVrGjaRXePdkWJyMCoLPLYkkQNe460seb11qijiEiOU1nksQ/PmsjIkjjLX9SuKBEZGJVFHhteEmfRhZN4csNu2jp1JVoR6T+VRZ6rm1/D0bYunt28L+ooIpLDVBZ5buF546koK9FZUSIyICqLPFcUj3H13Gp+/so+jrR1Rh1HRHKUyqIA1M+voaOrh6de1k2RRKR/VBYFYN7kcqaMH0GDrhUlIv2ksigAp65E+9/bWtlzuC3qOCKSg1QWBaI+UY07PL5eB7pFpO9UFgXi3Moy5k4u164oEekXlUUBqUvUsHHXEbbuOxp1FBHJMSqLAnLN3Cpihn5zISJ9prIoIBNGl/Le8ypYrpsiiUgfqSwKTF2imu0HTvDC9kNRRxGRHKKyKDCLL5rEsKKY7s8tIn2isigwo0qLuWLWRH7y0m46u3uijiMiOUJlUYDqEtUcON7Br17dH3UUEckRoZaFmS02s81mttXM7jzLPB8zs2Yz22hmD6eM32hmrwaPG8PMWWgunzmB8uHF2hUlIhkrCuuNzSwO3At8GNgBrDWzFe7enDLPdOCLwEJ3P2hmE4LxccBdQC3gwLpg2YNh5S0kJUUxrppTxfKmnZzo6GJESWh/DUQkT4S5ZbEA2Oru29y9A3gUqDttnluBe0+VgLufukPPImClux8IXlsJLA4xa8GpS1RzoqOblc17o44iIjmgz2VhZjEzG53BrDXAmynTO4KxVDOAGWb2azNbbWaL+7CsDMCCqeOoLi+l4UXtihKR9DIqCzN72MxGm9lIYAPQbGafH4TPLwKmA5cDS4H7zWxMpgub2W1m1mhmjS0tLYMQp3DEYsY1iWqee3U/rcfao44jIlku0y2L2e5+BKgHfgpMA25Is8xO4JyU6cnBWKodwAp373T314EtJMsjk2Vx9/vcvdbdaysrKzNcFTmlPlFDd4/zxMu7o44iIlku07IoNrNikmWxwt07SR547s1aYLqZTTOzEuA6YMVp8zSQ3KrAzCpI7pbaBjwNXGlmY81sLHBlMCaDaFbVaGZOHKVdUSKSVqZl8W3gN8BI4DkzmwIc6W0Bd+8Cbif5Jb8J+JG7bzSzu81sSTDb00CrmTUDq4DPu3urux8A7iFZOGuBu4MxGWR186t5YfshtreeiDqKiGQx6+8F5cysKCiErFBbW+uNjY1Rx8g5Ow6e4H1/s4o/v3IGt//29KjjiMgQM7N17l6bbr5MD3B/NjjAbWb2gJm9APz2gFNK5CaPHcElU8fSoCvRikgvMt0NdXNwgPtKYCzJg9tfDy2VDKm6RA1b9x1j465e9yyKSAHLtCws+PMq4HvuvjFlTHLcR+ZUURQzXf5DRM4q07JYZ2Y/I1kWT5vZKECXLM0TY0eWcPnMSlas30V3j3ZFici7ZVoWnwLuBC5x9xNACfDJ0FLJkKtL1LD3SDtrtrVGHUVEslBGZeHuPSR/GPdlM/t74L3u/lKoyWRIXTFrIiNL4ro/t4icUaZnQ30d+CzQHDzuMLO/DjOYDK3hJXEWXTSJJzfspq2zO+o4IpJlMt0NdRXwYXd/0N0fJHkF2KvDiyVRqEvUcLSti2c370s/s4gUlL5cdTb1An/lgx1EorfwvPFUlJXQ8KJ2RYnIO2V615uvAS+a2SqSp8x+gOQBb8kjRfEYV8+t5uE12zl8spPy4cVRRxKRLJHpAe5HgMuAfwd+DPwPd/9hmMEkGvXza+jo7uGpDboSrYi8rdeyMLP3nHoAVSQvKb4DqA7GJM/Mm1zO1PEjtCtKRN4h3W6of+jlNUfXh8o7ZkZdooZ/euZV9hxuY1J5adSRRCQL9FoW7v6hoQoi2aMuUc23fv4qj6/fxa0fODfqOCKSBTI6wG1mv3+G4cPAy+6u8yzzzLmVZcydXE5D006VhYgAfbvcxzLg+uBxP/AF4Ndmlu72qpKD6hI1bNx1hK37jkYdRUSyQKZlUQTMcvePuvtHgdkkj1lcSrI0JM9cM6+KmKED3SICZF4W57j73pTpfcHYAaBz8GNJ1CaMKmXh+RUsX79TN0USkYzL4lkz+4mZ3WhmNwIrgrGRwKHw4kmU6hI1vHngJC9sPxh1FBGJWKZl8Rng34BE8HgI+Iy7H9cZU/lr0YUTGVYU05VoRSTjX3A78CvgGeDnwHOufRN5b1RpMVfMnshPXtpNZ7fudSVSyDK9RPnHgOeBa4GPAWvM7Nowg0l2qJtXzYHjHfzq1f1RRxGRCGV6IcEvkbxL3j4AM6sE/hN4LKxgkh0unzmB8uHFNDTt5EMXTIg6johEJNNjFrHTfnzX2odlJYeVFMW4ak4VP9u4l+PtXVHHEZGIZPqF/5SZPW1mN5nZTcATwJPhxZJsUp+o5mRnNyub96afWUTyUqYHuD8P3AfMDR73ubt+jFcgLpk6juryUhqadkYdRUQikukxC9z9xyTvZSEFJhYzliRquP+X29h/rJ2KsmFRRxKRIZbufhZHzezIGR5HzezIUIWU6NXPr6a7x3nyZd0USaQQ9VoW7j7K3Uef4THK3UcPVUiJ3gWTRjNz4igaXtSuKJFCpDOaJGN186t5YfshtreeiDqKiAwxlYVkbMm8agCW60C3SMFRWUjGJo8dwYKp42ho0pVoRQqNykL6pG5+Na+1HGfjLp3fIFJIVBbSJx+ZU0Vx3HSgW6TAqCykT8aMKOGDMybw+Eu76O7RriiRQqGykD6rS1Sz90g7a7a1Rh1FRIZIqGVhZovNbLOZbTWzO8/w+k1m1mJmTcHjlpTXulPGV4SZU/rmilkTGVkS1+U/RApIaGVhZnHgXuB3gdnAUjObfYZZf+juieCxLGX8ZMr4krBySt8NL4mz6KJJ/PTlPbR1dkcdR0SGQJhbFguAre6+zd07gEeBuhA/T4ZQfaKGo+1drHplX/qZRSTnhVkWNcCbKdM7grHTfdTMXjKzx8zsnJTxUjNrNLPVZlYfYk7ph/eeN56KsmHaFSVSIKI+wP04MNXd5wIrgYdSXpvi7rXAJ4Bvmtl5py9sZrcFhdLY0tIyNIkFgKJ4jGvmVbHqlRYOn+iMOo6IhCzMstgJpG4pTA7G3uLure7eHkwuAy5OeW1n8Oc24Flg/ukf4O73uXutu9dWVlYObnpJqz5RQ0d3D09t1JVoRfJdmGWxFphuZtPMrAS4DnjHWU1mVpUyuQTYFIyPNbNhwfMKYCHQHGJW6Ye5k8uZOn4EDS/uijqKiIQstLJw9y7gduBpkiXwI3ffaGZ3m9mps5vuMLONZrYeuAO4KRifBTQG46uAr7u7yiLLmBl1iRpWv97KnsNtUccRkRBZvlwQrra21hsbG6OOUXBe33+cD/39s/zlVRdw2wfedVhJRLKcma0Ljg/3KuoD3JLjplWMZN7kcu2KEslzKgsZsLpEDc27j/Dq3qNRRxGRkKgsZMCunldFzNBvLkTymMpCBmzCqFIWnl/B8qZduimSSJ5SWcigqEvUsOPgSV7YfjDqKCISApWFDIpFF05kWFFMB7pF8pTKQgbFqNJirpg9kSde3k1nd0/UcURkkKksZNDUJ2o4cLyDX76q63SJ5BuVhQyaD86oZMyIYu2KEslDKgsZNCVFMa6aU8XK5r0cb++KOo6IDCKVhQyqunnVnOzsZmXz3qijiMggUlnIoLpk6jiqy0v1Az2RPKOykEEVixlLEjX88tX97D/Wnn4BEckJKgsZdPXzq+nucZ54STdFEskXKgsZdBdMGs0Fk0ZpV5RIHlFZSCjqEjW8uP0Qb7QejzqKiAwClYWEYkmiGoDlTfrNhUg+UFlIKGrGDGfB1HE0NO3UlWhF8oDKQkJTN7+abS3H2bjrSNRRRGSAVBYSmo/MqaI4bjS8qAPdIrlOZSGhGTOihA/OmMCK9bvo7tGuKJFcprKQUNXPr2bf0XZWb2uNOoqIDIDKQkJ1xayJlA0r0q4okRynspBQlRbHWXThJJ7asIe2zu6o44hIP6ksJHR1iWqOtnex6pV9UUcRkX5SWUjo3nveeCrKhunyHyI5TGUhoSuKx7hmXhWrXmnh8InOqOOISD+oLGRI1Cdq6Oju4acbdCVakVykspAhMXdyOdMqRmpXlEiOUlnIkDAz6hLVrHn9ALsPn4w6joj0kcpChkxdogZ3WKEr0YrkHJWFDJlpFSOZN7mcBpWFSM5RWciQqkvUsGn3EbbsPRp1FBHpA5WFDKmr51URM1iuA90iOUVlIUNqwqhSFp5fwfKmXbopkkgOUVnIkKtP1LDj4EnWvXEw6igikiGVhQy5RRdNorQ4pt9ciOSQUMvCzBab2WYz22pmd57h9ZvMrMXMmoLHLSmv3WhmrwaPG8PMKUOrbFgRV8yayBMv7aazuyfqOCKSgdDKwsziwL3A7wKzgaVmNvsMs/7Q3RPBY1mw7DjgLuBSYAFwl5mNDSurDL26RA0HT3Ty3JaWqKOISAbC3LJYAGx1923u3gE8CtRluOwiYKW7H3D3g8BKYHFIOSUCH5xRyZgRxSzXby5EckKYZVEDvJkyvSMYO91HzewlM3vMzM7py7JmdpuZNZpZY0uL/g81l5QUxbhqThUrm/dyvL0r6jgikkbUB7gfB6a6+1ySWw8P9WVhd7/P3WvdvbaysjKUgBKe+kQNJzu7+VnznqijiEgaYZbFTuCclOnJwdhb3L3V3duDyWXAxZkuK7mvdspYasYMp+FF7YoSyXZhlsVaYLqZTTOzEuA6YEXqDGZWlTK5BNgUPH8auNLMxgYHtq8MxiSPxGLGkkQ1v9q6n/3H2tMvICKRCa0s3L0LuJ3kl/wm4EfuvtHM7jazJcFsd5jZRjNbD9wB3BQsewC4h2ThrAXuDsYkz9QnaujucX6yXlsXItnM8uWSC7W1td7Y2Bh1DOmHxd98jtLiOA2fWRh1FJGCY2br3L023XxRH+AWoS5RQ9Obh3ij9XjUUUTkLFQWErkliWoA/eZCJIupLCRyNWOGs2DaOBqadupKtCJZSmUhWaE+UcO2luNs2Hkk6igicgYqC8kKV82ZRHHcdCVakSylspCsMGZECZfPnMDj63fR3aNdUSLZRmUhWaMuUc2+o+3892utUUcRkdOoLCRrXDFrImXDinR/bpEspLKQrFFaHGfRhZN4asMe2jq7o44jIilUFpJV6udXc7S9i2de2Rd1FBFJobKQrPLe8yqoHDWMhhe1K0okm6gsJKvEY8Y1c6t5dnMLh090Rh1HRAIqC8k6dYlqOrp7eHLD7qijiEigKOoAIqebO7mcaRUj+ceVW/j5pr2MKi1mVGlR8Eg+LxtWxOi3xoOx0iLKSoqIxSzqVRDJOyoLyTpmxucXzeQ7//Ubdh5q42jbUY62dXGsvSvtD/bMoKyk6F0lklo4o1MK5/Tx5FgRRXFtdIukUllIVrpqThVXzal6x5i7c7Kzm6NtXRxt6+RIW1eyRILp1PFj7W+PtR7r4Df7jydfb++io6sn7ecPL46/a2tmVGkRo4a9vTVT9lbJpM7zduGUFsfD+sfzLu6OOzjQ407PqWlPmQa8h7een3W+My3nTk8w9s7PPVOW06ZJv0za90zzHqcvcqYLUr57nt5z5JKRw+JcMGl0qJ+hspCcYWaMKCliREkRE0eX9vt92rtOFU6yUI61dQXF0/nW+LH2t58fCcZ3HToZlFAXJzrS/w6kJB57q2TiMQu+dM/8pfz2ePDFTzDd470vx9vLS+FKnDMm9JuHqSyk4AwrijOsLE5F2bB+v0dXd89bxXE0tWhSSiZ1vNsdA2JmxCxZfGYp0xixWDDe23xnWi7D+U5//azLnfr82OnLJZ+ffkTIzniIyHqd593vYb2+fsb3eNd7pvmQDD43V40qDf+rXGUh0g9F8RhjRpQwZkRJ1FFEhoSO4omISFoqCxERSUtlISIiaaksREQkLZWFiIikpbIQEZG0VBYiIpKWykJERNKyM11DJReZWQvwxgDeogLYP0hxopQv6wFal2yVL+uSL+sBA1uXKe5emW6mvCmLgTKzRnevjTrHQOXLeoDWJVvly7rky3rA0KyLdkOJiEhaKgsREUlLZfG2+6IOMEjyZT1A65Kt8mVd8mU9YAjWRccsREQkLW1ZiIhIWiqLgJndY2YvmVmTmf3MzKqjztRfZvZ3ZvZKsD7/YWZjos7UX2b2B2a20cx6zCznzlwxs8VmttnMtprZnVHnGQgze9DM9pnZhqizDISZnWNmq8ysOfi79dmoM/WXmZWa2fNmtj5Yl6+E9lnaDZVkZqPd/Ujw/A5gtrv/ScSx+sXMrgSecfcuM/sbAHf/QsSx+sXMZgE9wLeBP3f3xogjZczM4sAW4MPADmAtsNTdmyMN1k9m9gHgGPBdd78o6jz9ZWZVQJW7v2Bmo4B1QH0u/nux5K3+Rrr7MTMrBn4FfNbdVw/2Z2nLInCqKAIjeff93XOGu//M3buCydXA5CjzDIS7b3L3zVHn6KcFwFZ33+buHcCjQF3EmfrN3Z8DDkSdY6Dcfbe7vxA8PwpsAmqiTdU/nnQsmCwOHqF8d6ksUpjZX5nZm8D1wP+OOs8guRn4adQhClQN8GbK9A5y9EspX5nZVGA+sCbaJP1nZnEzawL2ASvdPZR1KaiyMLP/NLMNZ3jUAbj7l9z9HOAHwO3Rpu1dunUJ5vkS0EVyfbJWJusiMtjMrAz4MfC50/Ys5BR373b3BMk9CAvMLJRdhEVhvGm2cvcrMpz1B8CTwF0hxhmQdOtiZjcBVwO/41l+YKoP/15yzU7gnJTpycGYRCzYv/9j4Afu/u9R5xkM7n7IzFYBi4FBPwmhoLYsemNm01Mm64BXosoyUGa2GPgLYIm7n4g6TwFbC0w3s2lmVgJcB6yIOFPBCw4KPwBscvdvRJ1nIMys8tTZjmY2nOTJFKF8d+lsqICZ/RiYSfLMmzeAP3H3nPy/QDPbCgwDWoOh1Tl8ZtfvAf8MVAKHgCZ3XxRtqsyZ2VXAN4E48KC7/1XEkfrNzB4BLid5hdO9wF3u/kCkofrBzN4H/BJ4meR/7wB/6e5PRpeqf8xsLvAQyb9fMeBH7n53KJ+lshARkXS0G0pERNJSWYiISFoqCxERSUtlISIiaaksREQkLZWFSB+Y2bH0c/W6/GNmdm7wvMzMvm1mr5nZOjN71swuNbMSM3vOzArqR7OS3VQWIkPEzC4E4u6+LRhaRvLCfNPd/WLgk0BFcNHBnwMfjyapyLupLET6wZL+LriG1ctm9vFgPGZm/xLcT2SlmT1pZtcGi10PLA/mOw+4FPiyu/cAuPvr7v5EMG9DML9IVtBmrkj//D6QAOaR/EXzWjN7DlgITAVmAxNIXv76wWCZhcAjwfMLSf4avfss778BuCSU5CL9oC0Lkf55H/BIcMXPvcAvSH65vw/4f+7e4+57gFUpy1QBLZm8eVAiHcHNeUQip7IQGTongdLg+UZgXnA3vbMZBrSFnkokAyoLkf75JfDx4MYzlcAHgOeBXwMfDY5dTCR54b1TNgHnA7j7a0Aj8JXgKqiY2VQz+0jwfDyw3907h2qFRHqjshDpn/8AXgLWA88AfxHsdvoxyTviNQPfB14ADgfLPME7y+MWYCKw1cw2AN8hebczgA8F84tkBV11VmSQmVmZux8Ltg6eBxYsZAemAAAAbUlEQVS6+57gfgOrgumzHdg+9R7/Dtzp7luGILJIWjobSmTw/SS4IU0JcE+wxYG7nzSzu0jeh3v72RYObpTUoKKQbKItCxERSUvHLEREJC2VhYiIpKWyEBGRtFQWIiKSlspCRETSUlmIiEha/x9FE90wXwN9LwAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# scores_：dict with classes as the keys, and the values as the grid of scores obtained during cross-validating each fold,\n",
    "# Each dict value has shape (n_folds, len(Cs))\n",
    "Cs = [1e-3, 1e-2, 1e-1, 1, 10, 100, 1000]\n",
    "n_Cs = len(Cs)\n",
    "#Scores only 1 row to show the \n",
    "scores =  np.zeros((1,n_Cs))\n",
    "scores[0][:] = np.mean(lrcv_L1.scores_[1],axis = 0)\n",
    "    \n",
    "mse_mean = -np.mean(scores, axis = 0)\n",
    "plt.plot(np.log10(Cs), mse_mean.reshape(n_Cs,1)) \n",
    "#plt(np.log10(reg.Cs)*np.ones(3), [0.28, 0.29, 0.30])\n",
    "plt.xlabel('log(C)')\n",
    "plt.ylabel('logloss')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([0.69314718, 0.63596759, 0.47821638, 0.47497881, 0.47535088,\n",
       "       0.47540207, 0.47540823])"
      ]
     },
     "execution_count": 16,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "mse_mean"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "1 0.47497880655956004\n"
     ]
    }
   ],
   "source": [
    "best_C = np.argmin(mse_mean)\n",
    "best_score = np.min(mse_mean)\n",
    "print(Cs[best_C], best_score)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[ 0.40685268,  1.11908244, -0.09951999,  0.03399565, -0.06084834,\n",
       "         0.63683915,  0.29773472,  0.12255592,  0.12627302]])"
      ]
     },
     "execution_count": 18,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "lrcv_L1.coef_"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## L1正则 Solver='libliner' scoring = 'accuracy' (正确率)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "LogisticRegressionCV(Cs=[0.001, 0.01, 0.1, 1, 10, 100, 1000],\n",
       "           class_weight=None, cv=5, dual=False, fit_intercept=True,\n",
       "           intercept_scaling=1.0, max_iter=100, multi_class='ovr',\n",
       "           n_jobs=None, penalty='l1', random_state=None, refit=True,\n",
       "           scoring='accuracy', solver='liblinear', tol=0.0001, verbose=0)"
      ]
     },
     "execution_count": 19,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "from sklearn.linear_model import LogisticRegressionCV\n",
    "\n",
    "Cs = [1e-3, 1e-2, 1e-1, 1, 10, 100, 1000]\n",
    "#nCs = 1  #Cs values are chosen in a logarithmic scale between 1e-4 and 1e4.\n",
    "\n",
    "# LogisticRegressionCV比GridSearchCV快\n",
    "lrcv_L1_A = LogisticRegressionCV(Cs=Cs, cv = 5, scoring = 'accuracy', penalty='l1', solver='liblinear', multi_class='ovr')\n",
    "lrcv_L1_A.fit(X_train, y_train)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "{1: array([[0.64935065, 0.68181818, 0.75974026, 0.76623377, 0.76623377,\n",
       "         0.76623377, 0.76623377],\n",
       "        [0.64935065, 0.68181818, 0.74025974, 0.73376623, 0.73376623,\n",
       "         0.73376623, 0.73376623],\n",
       "        [0.64935065, 0.73376623, 0.77272727, 0.76623377, 0.76623377,\n",
       "         0.76623377, 0.76623377],\n",
       "        [0.65359477, 0.7254902 , 0.78431373, 0.79738562, 0.79738562,\n",
       "         0.79738562, 0.79738562],\n",
       "        [0.65359477, 0.7124183 , 0.75816993, 0.76470588, 0.76470588,\n",
       "         0.76470588, 0.76470588]])}"
      ]
     },
     "execution_count": 20,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "lrcv_L1_A.scores_"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYsAAAEKCAYAAADjDHn2AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzt3Xl4VfW97/H3l4SEIcyEQWYkQKmzEa0DFcp0O2hteyrWDnbQ0zpRrbTe2+e2vZ7b55zrgJUerQdb23qOlXqsemgPNQxCnVAJOFSGhBAZAsoQCBCGjN/7x15ptyGyV0JW1k7yeT3Pfp6stdfa+7MY9idr2Otn7o6IiMjJdIk7gIiIpD+VhYiIpKSyEBGRlFQWIiKSkspCRERSUlmIiEhKKgsREUlJZSEiIimpLEREJKXMuAO0loEDB/ro0aPjjiEi0q6sXbt2n7vnplquw5TF6NGjKSwsjDuGiEi7Ymbbwiynw1AiIpKSykJERFJSWYiISEoqCxERSUllISIiKaksREQkJZWFiIik1GG+ZyGSDo7X1LH3cBX7KqvYV1nNvsoqyiurqK6tjzuadGBD+nTnSxeOjPQ9VBYiJ+HuVFbV/u2Df1+jIvhgKVRTWVXb5OuYtXFw6VTOGdFXZSHS2urrnYPHathXWcXe4MO+vOGD//CJJVDVxF6BGfTrkcXAnCwG5mRz9vC+DAh+zs3JZmCvxM8Dc7Lp3zOLbl0zYthSkdajspAOobaunv1Hq//2YV9+5O8f/A2FsO9wYn55ZTW19X7Ca2R0MQb0DD7ke2Vzem4OA3tl/60Q/v7Ion/PLDIzdMpPOg+VhaStqto6yps43PPB3/4TH/77j1bjJ37+k5XZJfGbfk4WQ/t044xhvf/+od+oCPp270qXLjpeJNIUlYWkpflLi1jwfEmTz/XMygg+6LMZPaAn+aP7B4d//r5XMKBnFgN7ZdMrOxPTCQORU6aykLSzrfwID63awuUTcpn10SEMzMlmQE5WsIeQTfcsHf8XaWuRloWZzQYeADKAX7r7vzR6/n5gajDZAxjk7n2D50YCvwRGAA580t23RplX0sP9y4rJzDDu/vxZDOrdLe44IkKEZWFmGcCDwAygDFhjZovdfUPDMu5+W9LytwDnJr3EY8BP3X2ZmeUAulC9E9j0/iH+661d3DBlrIpCJI1EeTnHZKDE3UvdvRpYBFx5kuWvAZ4AMLNJQKa7LwNw90p3PxphVkkT9xYUk5OVyXc+fnrcUUQkSZRlMQzYkTRdFsw7gZmNAsYAzwezxgMVZva0mb1hZvcEeyrSga3bfoDlG3dzw5Sx9O2RFXccEUmSLheKzwGecve6YDoTuAy4A7gAGAtc13glM7vBzArNrHDv3r1tlVUi4O7c81wRA3pm8Y1Lx8QdR0QaibIsdpI4Od1geDCvKXMIDkEFyoA3g0NYtcCzwHmNV3L3he6e7+75ubkpxxuXNPZySTmrS8u5aeo4embrIj2RdBNlWawB8sxsjJllkSiExY0XMrOJQD9gdaN1+5pZQwNMAzY0Xlc6BnfnnoJNDOvbnWsvivb+NiLSMpGVRbBHcDNQAGwEnnT39WZ2l5ldkbToHGCR+9+/fxscjroDWGFmfwUMeCSqrBKvgvW7eavsIHOn55GdqVNTIunIvKl7JLRD+fn5XlhYGHcMaaa6emf2z16g3p2C707R/ZZE2piZrXX3/FTL6X+mxOrZN3ayeU8l35s5QUUhksb0v1NiU11bz/3LizljWG9mf3RI3HFE5CRUFhKbRWu2U3bgGPNmTdTdXkXSnMpCYnG0upYFK0qYPKY/U/IGxh1HRFJQWUgsfvPKVvZVVvH9WRN0C3GRdkBlIW3u4NEaHl61hWkTB5E/un/ccUQkBJWFtLmFL27h0PFa7pg5Ie4oIhKSykLa1J7Dx3n0pa185uzTmHRa77jjiEhIKgtpUw+t3EJ1XT23zxgfdxQRaQaVhbSZHfuP8vhr2/hi/nDGDOwZdxwRaQaVhbSZB1Zsxsy49RN5cUcRkWZSWUib2Lz7ME+vK+OrF41iaJ/ucccRkWZSWUibmL+smO5dM7hx6ri4o4hIC6gsJHJvl1Xw53fe51uXjaV/Tw2XKtIeqSwkcvcUFNGvR1e+dZmGSxVpr1QWEqnVW8p5cfM+brx8HL26dY07joi0kMpCItMwXOqQ3t34ysdGxR1HRE6BykIis2LjHtZtr+DWT+TRrauGSxVpz1QWEon6eufepUWMHtCDf8gfHnccETlFKguJxB/f3sWm9w9z24zxdNVwqSLtnv4XS6urqatn/rJiJg7pxWfOOi3uOCLSClQW0uqeLNzBtvKjzJs1QcOlinQQKgtpVcdr6liwYjPnj+rHtImD4o4jIq1EZSGt6rHVW9l9qIp5Gi5VpENRWUirOXy8hodWbWHK+FwuGjsg7jgi0opUFtJqHnnxXSqO1jBPw6WKdDiRloWZzTazIjMrMbM7m3j+fjN7M3gUm1lFo+d7m1mZmf1rlDnl1JVXVvGrF0v55JlDOHN4n7jjiEgry4zqhc0sA3gQmAGUAWvMbLG7b2hYxt1vS1r+FuDcRi/zT8ALUWWU1vPQqi0cq6nTcKkiHVSUexaTgRJ3L3X3amARcOVJlr8GeKJhwszOBwYDSyPMKK1gV8Ux/v3VbXz+vOGMG9Qr7jgiEoEoy2IYsCNpuiyYdwIzGwWMAZ4PprsA9wF3RJhPWsmCFZvBYe50DZcq0lGlywnuOcBT7l4XTN8ILHH3spOtZGY3mFmhmRXu3bs38pByotK9lfzn2jK+dOFIhvfrEXccEYlIZOcsgJ3AiKTp4cG8pswBbkqa/hhwmZndCOQAWWZW6e4fOEnu7guBhQD5+fneWsElvPnLisnO7MJNGi5VpEOLsizWAHlmNoZEScwBvtR4ITObCPQDVjfMc/drk56/DshvXBQSv/W7DvKnt9/j5qnjyO2VHXccEYlQZIeh3L0WuBkoADYCT7r7ejO7y8yuSFp0DrDI3bVn0M7cW1BEn+5duX7K2LijiEjEotyzwN2XAEsazftRo+mfpHiN3wC/aeVocorWbN3PyqK9/GD2RPp013CpIh1dupzglnbE3bnnuSJye2Vz3cWj444jIm1AZSHNtqp4L69v3c+t08bRPUvDpYp0BioLaZb6eufegiJG9O/O1ReMjDuOiLQRlYU0y5J33mP9rkPcNn08WZn65yPSWeh/u4RWW1fP/KXFjB+cw5XnNPllfBHpoFQWEtof1pVRuu8I35s5gQwNlyrSqagsJJTjNXU8sHwzZ4/oy8xJg+OOIyJtTGUhoTz+2nZ2HTzO9zVcqkinpLKQlCqranloZQmXjBvAJeMGxh1HRGKgspCUHn3pXcqPVHOHhksV6bRUFnJSB45U88gLpcycNJhzR/aLO46IxERlISf18F+2UFldy/e0VyHSqaks5EPtPnSc37yylavOGcaEIRouVaQzS1kWZtbDzP63mT0STOeZ2aejjyZxW7BiM3X1znenj487iojELMyexa+BKhKj10FiIKP/G1kiSQvbyo/w+zU7uGbySEYO0HCpIp1dmLI43d3vBmoA3P0ooAvtO7j7lxWTmWHcMk3DpYpIuLKoNrPugAOY2ekk9jSkg9r0/iH+661dXHfxGAb17hZ3HBFJA2FGyvsx8BwwwsweBy4BrosylMTr3oJicrIz+fbHNVyqiCSctCwscV+HTcDngItIHH6a6+772iCbxGDd9gMs37ibO2aOp2+PrLjjiEiaOGlZuLub2RJ3PxP47zbKJDFpGC51YE4WX79kTNxxRCSNhDlnsc7MLog8icTu5ZJyVpeWc9PUcfTMDnOEUkQ6izCfCBcC15rZNuAIiUNR7u5nRZpM2pS7c0/BJob17c6XLtRwqSLyQWHKYlbkKSR2Bet381bZQe7+wllkZ2bEHUdE0kzKw1Duvg3oC3wmePQN5kkHUVfv3Le0iNNze/K5czVcqoicKMztPuYCjwODgsd/mNktUQeTtvPsGzvZvKeS782cQGaGbhcmIicK88nwTeBCd/+Ru/+IxCW014d5cTObbWZFZlZiZnc28fz9ZvZm8Cg2s4pg/jlmttrM1pvZ22Z2dXM2SsKrrq3n/uXFnDGsN7M/OiTuOCKSpsKcszCgLmm6jhC3+zCzDOBBYAZQBqwxs8XuvqFhGXe/LWn5W4Bzg8mjwFfdfbOZnQasNbMCd68IkVeaYdGa7ZQdOMZPrzqTLl10FxcRaVqYsvg18JqZPRNMfxb4VYj1JgMl7l4KYGaLgCuBDR+y/DUkvi2Ouxc3zHT3XWa2B8gFVBat6Gh1LQtWlDB5TH+m5Gm4VBH5cCnLwt3nm9kq4NJg1tfd/Y0Qrz0M2JE0XUbiMtwTmNkoYAzwfBPPTQaygC0h3lOa4TevbGVfZRUPf/k8El/WFxFpWsqyMLOLgPXuvi6Y7m1mF7r7a62YYw7wlLsnH+7CzIYC/w58zd3rm8h2A3ADwMiR+m5Acxw8WsPDq7YwbeIg8kf3jzuOiKS5MCe4fwFUJk1XBvNS2QmMSJoeHsxryhzgieQZZtabxC1Gfujurza1krsvdPd8d8/Pzc0NEUkaLHxxC4eO13KHhksVkRDClIW5uzdMBL/hhznXsQbIM7MxZpZFohAWn/DiZhOBfsDqpHlZwDPAY+7+VIj3kmbYc/g4j760lc+cfRqTTusddxwRaQfClEWpmd1qZl2Dx1ygNNVK7l4L3AwUABuBJ919vZndZWZXJC06B1iUXEjAF4EpwHVJl9aeE3qr5KQeWrmF6rp6bp+h4VJFJBz74Gd0EwuYDQIWANNIDIC0Aviuu++JPl54+fn5XlhYGHeMtLdj/1Gm3beKL5w/nH/+nG7vJdLZmdlad89PtVyYq6H2kPjtXzqAB1Zsxsy49RN5cUcRkXYkzO0+7g6ugOpqZivMbK+Zfbktwknr2rz7ME+vK+OrF41iaJ/ucccRkXYkzDmLme5+CPg0sBUYB8yLMpREY/6yYnpkZXLj1HFxRxGRdiZMWTQcqvoU8J/ufjDCPBKRt8sq+PM77/Oty8bQv6eGSxWR5glzCeyfzGwTcAz4jpnlAsejjSWt7Z6CIvr16Mo3L9VwqSLSfGHGs7gTuBjId/caEjf5uzLqYNJ6Vm8p58XN+7hp6jh6desadxwRaYdCDbTs7vuTfj5CYnhVaQcahksd0rsbX75oVNxxRKSd0kg3HdyKjXtYt72CudPz6NZVw6WKSMuoLDqw+nrn3qVFjB7Qgy+cPzzuOCLSjoX5nsXTZvYpM1OxtDN/fHsXm94/zO0zJ9BVw6WKyCkI8wnyEPAlYLOZ/YuZ6Tal7UBNXT3zlxXzkaG9+fSZQ+OOIyLtXJiroZa7+7XAeSS+lLfczF4xs6+bmS6tSVNPFu5gW/lR5s0ar+FSReSUhTo2YWYDgOuAbwFvAA+QKI9lkSWTFjteU8eCFZs5f1Q/pk4YFHccEekAwoyU9wwwgcSIdZ9x9/eCp35vZrrNaxp6bPVWdh+qYsGcczVcqoi0ijDfs1jg7iubeiLMbW2lbR0+XsNDq7YwZXwuF44dEHccEekgwhyGmmRmfRsmzKyfmd0YYSY5BY+8+C4VR2uYp+FSRaQVhSmL6929omHC3Q8A10cXSVqqvLKKX71YyifPHMKZw/vEHUdEOpAwZZFhSQe+zSwD0G1L09BDq7ZwrKaO22dor0JEWleYcxbPkTiZ/W/B9D8G8ySN7Ko4xr+/uo3PnzeccYNy4o4jIh1MmLL4AYmC+E4wvQz4ZWSJpEUWrNgMDnOna7hUEWl9Ycbgrgd+ETwkDZXureQ/15bxlYtGMbxfj7jjiEgHFOZ7FnnAPwOTgG4N8919bIS5pBnuX76Z7Mwu3KThUkUkImFOcP+axF5FLTAVeAz4jyhDSXjrdx3kj2/t4huXjCG3V3bccUSkgwpTFt3dfQVg7r7N3X9CYjxuSQP3LS2mT/euXD9FO3oiEp0wJ7irgtuTbzazm4GdgC63SQOFW/fz/KY9/GD2RPp01z0dRSQ6YfYs5gI9gFuB84EvA1+LMpSk5u7c/VwRub2yue7i0XHHEZEO7qRlEXwB72p3r3T3Mnf/urt/3t1fDfPiZjbbzIrMrMTM7mzi+fvN7M3gUWxmFUnPfc3MNgcPlVMjL2zex+tb93PrtHF0z9JwqSISrZMehnL3OjO7tCUvHBTNg8AMoAxYY2aL3X1D0uvflrT8LcC5wc/9gR8D+YADa4N1D7QkS0dTX+/cU7CJ4f26c/UFI+OOIyKdQJjDUG+Y2WIz+4qZfa7hEWK9yUCJu5e6ezWwCLjyJMtfAzwR/DwLWObu+4OCWAbMDvGencJz69/nnZ2HuG36eLIyNVyqiEQvzAnubkA5MC1pngNPp1hvGLAjaboMuLCpBc1sFDAGeP4k6w5rYr0bgBsARo7sHL9h19bVc+/SIvIG5fDZc0/4IxERiUSYb3B/vQ1yzAGecve65qzk7guBhQD5+fkeRbB08/QbOynde4SHv3w+GRouVUTaSJhvcP+axJ7EB7j7N1KsuhMYkTQ9PJjXlDnATY3WvbzRuqtSvF+HV1VbxwPLN3P28D7M+ujguOOISCcS5oD3n4D/Dh4rgN5AZYj11gB5ZjbGzLJIFMLixguZ2USgH7A6aXYBMDMYaKkfMDOY16n97rXt7Kw4xrxZEzVcqoi0qTCHof6QPG1mTwAvhVivNvgSXwGQATzq7uvN7C6g0N0bimMOsMjdPWnd/Wb2TyQKB+Aud98faos6qCNVtTy4soSPjR3AJeM0XKqItK0wJ7gbywMGhVnQ3ZcASxrN+1Gj6Z98yLqPAo+2IF+H9OuX32VfZTULvzpBexUi0ubCnLM4zAfPWbxPYowLaSMVR6v5txdKmf6RwZw3sl/ccUSkEwpzGKpXWwSRD/fwX0qprKrljlnj444iIp1UyhPcZnaVmfVJmu5rZp+NNpY02HPoOL955V2uPPs0Jg7pHXccEemkwlwN9WN3P9gw4e4VJG7FIW3g58+XUFvn3DZDexUiEp8wZdHUMi05MS7NtL38KE+8vp2rLxjBqAE9444jIp1YmLIoNLP5ZnZ68JgPrI06mMDPVhST0cW4ZVpe3FFEpJMLUxa3ANXA70ncDPA4H/y2tUSgePdhnnljJ9ddPJohfbqlXkFEJEJhroY6ApwwFoVE676lReRkZfLtj58edxQRkVBXQy0zs75J0/3MrNPfeiNKb+6ooGD9bq6fMpZ+PbPijiMiEuow1MDgCigAgvElQn2DW1rm3oIiBvTM4huXjok7iogIEK4s6s3sb4NFBGNPdIrbgcfhlZJ9vFSyjxunjiMnWxediUh6CPNp9EPgJTP7C2DAZQQDDknrcnfuLihiaJ9uXHth5xjMSUTahzAnuJ8zs/OAi4JZ33X3fdHG6pyWb9zDmzsq+JfPnUm3rhlxxxER+ZuwxznqgD0khlidZGa4+wvRxep86uqdewuKGDOwJ184f3jccUREPiDMXWe/BcwlMVrdmyT2MFbzwTG55RT98a1dFO0+zM+vOZfMjDCnkkRE2k6YT6W5wAXANnefCpwLVJx8FWmO6tp65i8rZtLQ3nzqzKFxxxEROUGYsjju7scBzCzb3TcBE6KN1bk8WbiD7fuPMm/WBLp00cBGIpJ+wpyzKAu+lPcssMzMDgDboo3VeRyrrmPBis3kj+rH5RNy444jItKkMFdDXRX8+BMzWwn0AZ6LNFUn8tjqrew5XMW/fuk8DZcqImmrWd/6cve/RBWkMzp0vIZf/GULHx+fy+Qx/eOOIyLyoXTZTYx++UIpFUdrmDdLp4BEJL2pLGKyr7KKX770Lp86cyhnDOuTegURkRipLGLy0MotHK+p4/aZGi5VRNKfyiIGOyuO8R+vbuML5w/n9NycuOOIiKSksojBguWbAZg7XXsVItI+RFoWZjbbzIrMrMTMmhxtz8y+aGYbzGy9mf0uaf7dwbyNZrbAOsh1pVv2VvLUujKuvWgkw/p2jzuOiEgokQ2YYGYZwIPADKAMWGNmi919Q9IyecD/BC5x9wNmNiiYfzFwCXBWsOhLwMeBVVHlbSv3LysmO7MLN00dF3cUEZHQotyzmAyUuHupu1cDi4ArGy1zPfBgMPoe7r4nmO8k7nCbBWQDXYHdEWZtE+/sPMif3n6Pb146hoE52XHHEREJLcqyGAbsSJouC+YlGw+MN7OXzexVM5sN4O6rgZXAe8GjwN03Rpi1Tdy3tIg+3bvyrcvGxh1FRKRZ4h63MxPIAy4ncQv0F8zsTGAg8JFgHiTuSXWZu7+YvLKZ3UAwat/Ikek9styarftZWbSXO//HRPp07xp3HBGRZolyz2InMCJpengwL1kZsNjda9z9XaCYRHlcBbzq7pXuXgn8GfhY4zdw94Xunu/u+bm56XsTPnfn7uc2kdsrm699bHTccUREmi3KslgD5JnZGDPLAuYAixst8yyJvQrMbCCJw1KlwHbg42aWaWZdSZzcbreHoVYV72XN1gPcOm0c3bM0XKqItD+RlYW71wI3AwUkPuifdPf1ZnaXmV0RLFYAlJvZBhLnKOa5eznwFLAF+CvwFvCWu/8xqqxRqg+GSx3RvztXX5Deh8pERD5MpOcs3H0JsKTRvB8l/ezA7cEjeZk64B+jzNZW/vzO+6zfdYj5XzybrEx9B1JE2id9ekWotq6e+5YVMX5wDlee0/hCMBGR9kNlEaGn1+2kdO8RvjdzAhkaLlVE2jGVRUSO19Txs+XFnD2iLzMnDY47jojIKVFZROR3r21n18HjfH/WBA2XKiLtnsoiAkeqanlwZQkXnz6AS8YNjDuOiMgpU1lE4NGX3qX8SLWGSxWRDkNl0coqjlaz8IVSZkwazLkj+8UdR0SkVagsWtkv/rKFyupa7pipvQoR6ThUFq1o96Hj/PaVrXz2nGFMGNIr7jgiIq1GZdGKfv78ZmrrnO9Oz4s7iohIq1JZtJLt5UdZ9PoO5kwewagBPeOOIyLSqlQWreRny4vJzDBumaa9ChHpeFQWraDo/cM88+ZOvnbxaAb37hZ3HBGRVqeyaAX3LS0iJyuTb085Pe4oIiKRUFmcoje2H2Dpht3cMGUs/XpmxR1HRCQSKotTdO/SIgb0zOLrl46JO4qISGRUFqfg5ZJ9vFxSzo1Tx5GTHek4UiIisVJZtJC7c3dBEaf16ca1F2q4VBHp2FQWLbRsw27e2lHB3Ol5dOuaEXccEZFIqSxaoK7euXdpEWMH9uTz5w2PO46ISORUFi2w+K2dFO+u5PaZ48nM0B+hiHR8+qRrpuraeuYvK2bS0N588oyhcccREWkTKotm+n3hDnbsP8a82RPo0kXDpYpI56CyaIZj1XX8fMVmLhjdj8vH58YdR0SkzagsmuG3q7ey53AV82ZNxEx7FSLSeagsQjp0vIZfrNrC5RNymTymf9xxRETaVKRlYWazzazIzErM7M4PWeaLZrbBzNab2e+S5o80s6VmtjF4fnSUWVN55IVSDh6r0XCpItIpRXaPCjPLAB4EZgBlwBozW+zuG5KWyQP+J3CJux8ws0FJL/EY8FN3X2ZmOUB9VFlT2VdZxa9eepdPnTWUM4b1iSuGiEhsotyzmAyUuHupu1cDi4ArGy1zPfCgux8AcPc9AGY2Cch092XB/Ep3Pxph1pN6cGUJVbX13D5jfFwRRERiFWVZDAN2JE2XBfOSjQfGm9nLZvaqmc1Oml9hZk+b2Rtmdk+wp/IBZnaDmRWaWeHevXsj2YidFcd4/NXtfOG84ZyemxPJe4iIpLu4T3BnAnnA5cA1wCNm1jeYfxlwB3ABMBa4rvHK7r7Q3fPdPT83N5pLWR9YXgzArdM1XKqIdF5RlsVOYETS9PBgXrIyYLG717j7u0AxifIoA94MDmHVAs8C50WYtUlb9lby1NoyvnzRKIb17d7Wby8ikjaiLIs1QJ6ZjTGzLGAOsLjRMs+S2KvAzAaSOPxUGqzb18wadhemARtoY/OXFdOtawY3TtVwqSLSuUVWFsEewc1AAbAReNLd15vZXWZ2RbBYAVBuZhuAlcA8dy939zoSh6BWmNlfAQMeiSprU97ZeZD/fvs9vnnpGAbmZLflW4uIpB1z97gztIr8/HwvLCxstde77tev88b2Cl78wVR6d+vaaq8rIpJOzGytu+enWi7uE9xp6fV397OqaC/fufx0FYWICCqLE7g79xRsYlCvbL72sdFxxxERSQsqi0ZWFe1lzdYD3PKJPLpnabhUERFQWXxAfb1zT0ERI/v34Or8EalXEBHpJFQWSZa88x4b3jvEbTPyyMrUH42ISAN9IgZq6+qZv7SY8YNzuOLsxnclERHp3FQWgT+sK6N03xHumDmBDA2XKiLyASoL4HhNHT9bvplzRvRlxqTBcccREUk7Kgvg8de2897B43x/1gQNlyoi0oROXxaVVbU8uLKES8YN4OJxA+OOIyKSliIbKa+9OFpVy+TR/fn25bpZoIjIh+n0ZTGodzce/sr5cccQEUlrnf4wlIiIpKayEBGRlFQWIiKSkspCRERSUlmIiEhKKgsREUlJZSEiIimpLEREJCVz97gztAoz2wtsO4WXGAjsa6U4ceoo2wHalnTVUbalo2wHnNq2jHL33FQLdZiyOFVmVuju+XHnOFUdZTtA25KuOsq2dJTtgLbZFh2GEhGRlFQWIiKSksri7xbGHaCVdJTtAG1Luuoo29JRtgPaYFt0zkJERFLSnoWIiKSksgiY2T+Z2dtm9qaZLTWz0+LO1FJmdo+ZbQq25xkz6xt3ppYys38ws/VmVm9m7e7KFTObbWZFZlZiZnfGnedUmNmjZrbHzN6JO8upMLMRZrbSzDYE/7bmxp2ppcysm5m9bmZvBdvyfyJ7Lx2GSjCz3u5+KPj5VmCSu3875lgtYmYzgefdvdbM/h+Au/8g5lgtYmYfAeqBfwPucPfCmCOFZmYZQDEwAygD1gDXuPuGWIO1kJlNASqBx9z9jLjztJSZDQWGuvs6M+sFrAU+2x7/XszMgJ7uXmlmXYGXgLnu/mprv5f2LAINRRHoCbTbFnX3pe5eG0y+CgyPM8+pcPeN7l4Ud44WmgyUuHupu1cDi4AcRmpvAAACo0lEQVQrY87UYu7+ArA/7hynyt3fc/d1wc+HgY3AsHhTtYwnVAaTXYNHJJ9dKoskZvZTM9sBXAv8KO48reQbwJ/jDtFJDQN2JE2X0U4/lDoqMxsNnAu8Fm+SljOzDDN7E9gDLHP3SLalU5WFmS03s3eaeFwJ4O4/dPcRwOPAzfGmPblU2xIs80OglsT2pK0w2yLS2swsB/gD8N1GRxbaFXevc/dzSBxBmGxmkRwizIziRdOVu08PuejjwBLgxxHGOSWptsXMrgM+DXzC0/zEVDP+XtqbncCIpOnhwTyJWXB8/w/A4+7+dNx5WoO7V5jZSmA20OoXIXSqPYuTMbO8pMkrgU1xZTlVZjYb+D5whbsfjTtPJ7YGyDOzMWaWBcwBFsecqdMLTgr/Ctjo7vPjznMqzCy34WpHM+tO4mKKSD67dDVUwMz+AEwgceXNNuDb7t4ufws0sxIgGygPZr3ajq/sugr4OZALVABvuvuseFOFZ2afBH4GZACPuvtPY47UYmb2BHA5iTuc7gZ+7O6/ijVUC5jZpcCLwF9J/H8H+F/uviS+VC1jZmcBvyXx76sL8KS73xXJe6ksREQkFR2GEhGRlFQWIiKSkspCRERSUlmIiEhKKgsREUlJZSESITMbYmaLzGyLma01syVmNj7uXCLN1am+wS3SloIvfz0D/Nbd5wTzzgYGk7gbrUi7obIQic5UoMbdH26Y4e5vxZhHpMV0GEokOmeQGCtBpN1TWYiISEoqC5HorAfOjzuESGtQWYhE53kg28xuaJhhZmeZ2WUxZhJpEZWFSESCcUSuAqYHl86uB/4ZeD/eZCLNp7vOiohIStqzEBGRlFQWIiKSkspCRERSUlmIiEhKKgsREUlJZSEiIimpLEREJCWVhYiIpPT/AWj74NAY1600AAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# scores_：dict with classes as the keys, and the values as the grid of scores obtained during cross-validating each fold,\n",
    "# Each dict value has shape (n_folds, len(Cs))\n",
    "Cs = [1e-3, 1e-2, 1e-1, 1, 10, 100, 1000]\n",
    "n_Cs = len(Cs)\n",
    "#Scores only 1 row to show the \n",
    "scores =  np.zeros((1,n_Cs))\n",
    "scores[0][:] = np.mean(lrcv_L1_A.scores_[1],axis = 0)\n",
    "    \n",
    "mse_mean = np.mean(scores, axis = 0)\n",
    "plt.plot(np.log10(Cs), mse_mean.reshape(n_Cs,1)) \n",
    "#plt(np.log10(reg.Cs)*np.ones(3), [0.28, 0.29, 0.30])\n",
    "plt.xlabel('C')\n",
    "plt.ylabel('accuracy score')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([0.6510483 , 0.70706222, 0.76304219, 0.76566505, 0.76566505,\n",
       "       0.76566505, 0.76566505])"
      ]
     },
     "execution_count": 22,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "mse_mean"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "best_C = np.argmax(mse_mean)\n",
    "best_score = np.max(mse_mean)\n",
    "print(Cs[best_C], best_score)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[ 0.40689182,  1.11916565, -0.0995387 ,  0.03394157, -0.06095842,\n",
       "         0.636971  ,  0.29774113,  0.12257256,  0.12625301]])"
      ]
     },
     "execution_count": 23,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "lrcv_L1_A.coef_"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## L1正则 Solver='libliner' scoring = 'roc_auc' "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "LogisticRegressionCV(Cs=[0.001, 0.01, 0.1, 1, 10, 100, 1000],\n",
       "           class_weight=None, cv=5, dual=False, fit_intercept=True,\n",
       "           intercept_scaling=1.0, max_iter=100, multi_class='ovr',\n",
       "           n_jobs=None, penalty='l1', random_state=None, refit=True,\n",
       "           scoring='roc_auc', solver='liblinear', tol=0.0001, verbose=0)"
      ]
     },
     "execution_count": 24,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "from sklearn.linear_model import LogisticRegressionCV\n",
    "\n",
    "Cs = [1e-3, 1e-2, 1e-1, 1, 10, 100, 1000]\n",
    "#nCs = 1  #Cs values are chosen in a logarithmic scale between 1e-4 and 1e4.\n",
    "\n",
    "# LogisticRegressionCV比GridSearchCV快\n",
    "lrcv_L1_roc = LogisticRegressionCV(Cs=Cs, cv = 5, scoring = 'roc_auc', penalty='l1', solver='liblinear', multi_class='ovr')\n",
    "lrcv_L1_roc.fit(X_train, y_train)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "{1: array([[0.5       , 0.7562963 , 0.81462963, 0.81351852, 0.81425926,\n",
       "         0.81425926, 0.81407407],\n",
       "        [0.5       , 0.74898148, 0.79407407, 0.80111111, 0.80222222,\n",
       "         0.80203704, 0.80203704],\n",
       "        [0.5       , 0.82425926, 0.85222222, 0.84759259, 0.84555556,\n",
       "         0.84537037, 0.84537037],\n",
       "        [0.5       , 0.8240566 , 0.88471698, 0.88660377, 0.88566038,\n",
       "         0.8854717 , 0.8854717 ],\n",
       "        [0.5       , 0.80066038, 0.84075472, 0.84735849, 0.84830189,\n",
       "         0.84830189, 0.84830189]])}"
      ]
     },
     "execution_count": 25,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "lrcv_L1_roc.scores_"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYsAAAEKCAYAAADjDHn2AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzt3XuYXXV97/H3Z2YySQiX3CaAuYO5oSjgCFpEgQTI01aQ0+e08fIUrMqxpyDa1hZPe9SDp0/Vc7TaPtRKNWpbFanWnrRNuWQCcqeZUBSSPSEXwEyA2ZP7PXP7nj/2CuwMk+w1k1mzZ+/9eT3Pftjrt39r7e8KsL7Z6/f7rp8iAjMzsxOpK3cAZmY2+jlZmJlZSU4WZmZWkpOFmZmV5GRhZmYlOVmYmVlJThZmZlaSk4WZmZXkZGFmZiU1ZHlwSUuBrwP1wLci4ov9Pp8FfA+YmPS5LSJWSpoD5IANSdcnIuLjJ/quqVOnxpw5c4Y1fjOzard27drtEdFUql9myUJSPXAHcBXQDqyRtCIi1hd1+1Pg7oj4hqTzgJXAnOSzzRFxQdrvmzNnDq2trcMTvJlZjZD0Ypp+Wd6GuhjYFBFbIqILuAu4rl+fAE5P3p8BvJRhPGZmNkRZJovpwNai7fakrdjngQ9Jaqfwq+KWos/mSvpPST+TdFmGcZqZWQnlHuB+P/DdiJgB/Crw95LqgJeBWRFxIfD7wA8knd5/Z0k3SWqV1NrZ2TmigZuZ1ZIsk8U2YGbR9oykrdhHgLsBIuJxYBwwNSKORMSOpH0tsBmY3/8LIuLOiGiOiOamppLjM2ZmNkRZJos1wDxJcyU1AsuAFf36/BJYDCBpEYVk0SmpKRkgR9I5wDxgS4axmpnZCWQ2GyoieiTdDNxLYVrs8ohYJ+l2oDUiVgB/APytpE9RGOy+MSJC0ruB2yV1A33AxyNiZ1axmpnZialaVsprbm4OT501MxscSWsjorlUv0yL8swqSV9f0NXbR3dvHz29QXdvX7Id9BS97076dPcG3T3J+76i9/36dfUGpPlLmXTij1OcQ4lDoBRHKX2MgfvqODse06do72PbB9f/2OMX9TnucY7TP80fKun+7Ad1wEEcM80hp0wYy9I3n5X6u4fCycJGpee3H2DNCztff+HuCXr6Xnvf3dtX2O6JAS/URy/WPUWfdfUU9jl6se/q7aOnL+jty+5Xdqn/4avkB76VyQUzJzpZWG15ec8hvr5qI/+4tv24F28JGuvraKyvo6FejKmvY0x9HY0NdTTUJdsNdTTWi4a6OsY3Ft6Pqa+job6OMfWiMdmnoej9sdtiTEMdY+rqGNPw2neMqS9+f+x2Y32hb0Pda+/H1BdiOt7fuodTqVvKaRJSqS7F3xHHtBcfIwZsP14sx+t/7PFLfy8nccwTSZvHB5PwI+1RU3ZrqM++CsLJwkaF3Qe7+MaDm/nuYy8QAb/9ztl86B2zOXVsw+suyvV12V94K1GphDQ8+cp/9rXKycLK6lBXL8sffZ6/+dlm9h/p4foLp/OpJfOZOfmUcodmZkWcLKwsunv7+NGarfxly0by+46wZNE0Pn3NQhacdVq5QzOzAThZ2Ijq6wv+7ZmX+cp9G3hhx0HePmcSf/3Bi2ieM7ncoZnZCThZ2IiICB7euJ0v39vGs9v2suDM0/j2Dc1cuXDaiAz+mtnJcbKwzP18626+dE8bj23ewfSJ4/nqb76V6y6Y7oFqswriZGGZ2dy5n/977wb+/dlXmDKhkc+99zw+cMksxjbUlzs0MxskJwsbdsW1EuMa6vjkknl89LJzOHWs/3Mzq1T+v9eGTXGtRF8Ev/3O2fzeFW9k6qljyx2amZ0kJws7aa6VMKt+ThY2ZAPVSvzhNQtYeNbrFjU0swrnZGGD1r9Wonn2JO744EW83bUSZlXLycIG5eGNnXzpHtdKmNUaJwtL5edbd/Ple9t4dJNrJcxqkZOFndDmzv185b4NrHzGtRJmtczJwgb0yp7DfL3lOe5uda2EmTlZWD+7D3bxjZ9t5ruPulbCzF6TabKQtBT4OlAPfCsivtjv81nA94CJSZ/bImJl8tlngI8AvcAnIuLeLGOtdYe6evnOY8/zNw9uZp9rJcysn8yShaR64A7gKqAdWCNpRUSsL+r2p8DdEfENSecBK4E5yftlwJuANwCrJM2PiN6s4q1V3b193N26la+vcq2EmR1flr8sLgY2RcQWAEl3AdcBxckigKNXpTOAl5L31wF3RcQR4HlJm5LjPZ5hvDWlry9Y+ezLfOW+53h++wHXSpjZCWWZLKYDW4u224FL+vX5PHCfpFuACcCSon2f6Lfv9GzCrD0Pb+zky/ds4Jlte1wrYWaplHuA+/3AdyPiK5LeCfy9pDen3VnSTcBNALNmzcooxOrhWgkzG6osk8U2YGbR9oykrdhHgKUAEfG4pHHA1JT7EhF3AncCNDc3x7BFXmVcK2FmJyvLZLEGmCdpLoUL/TLgA/36/BJYDHxX0iJgHNAJrAB+IOmrFAa45wH/kWGsVcm1EmY2XDK7akREj6SbgXspTItdHhHrJN0OtEbECuAPgL+V9CkKg903RkQA6yTdTWEwvAf4Pc+ESs+1EmY23FS4Nle+5ubmaG1tLXcYZdXT28edD29xrYSZpSZpbUQ0l+rn+xFV5O7Wdr58zwYWL5zGp5e6VsLMho+TRRW5b/0rzJ5yCt+6odnTYM1sWNWVOwAbHge7enhs8w4WLzzTicLMhp2TRZV4ZON2unr6WLxoWrlDMbMq5GRRJVpyeU4b2+DHdZhZJpwsqkBfX7B6Q553L2iiscH/Ss1s+PnKUgWe2baHzn1HWLzQt6DMLBtOFlWgpS1PneDyBU4WZpYNJ4sq0JLr4KJZk5g8obHcoZhZlXKyqHAv7znEupf2snjRmeUOxcyqmJNFhVvdlgdgiafMmlmGnCwqXEsuz8zJ43njtFPLHYqZVTEniwp2qKuXRzdtd9W2mWXOyaKCPbppO0d6+lji8Qozy5iTRQVraevg1LENXDzXVdtmli0niwoVEbTk8rx7/lRXbZtZ5nyVqVDPbttLft8RrlzoW1Bmlj0niwq1KteBBFcsaCp3KGZWA5wsKtTqtjwXzZrEFK+rbWYjwMmiAnXsPcwz2/ZwpR8caGYjxMmiArXkjlZte7zCzEZGpslC0lJJGyRtknTbAJ//haSnk9dzknYXfdZb9NmKLOOsNKvbOpgxaTzzz3TVtpmNjIasDiypHrgDuApoB9ZIWhER64/2iYhPFfW/Bbiw6BCHIuKCrOKrVIe7e3lk03Z+q3mmq7bNbMRk+cviYmBTRGyJiC7gLuC6E/R/P/DDDOOpCo9t3s7h7j4/ZdbMRlSWyWI6sLVouz1pex1Js4G5wOqi5nGSWiU9Iel9x9nvpqRPa2dn53DFPaqtyuWZ0FjPJee4atvMRs5oGeBeBvw4InqL2mZHRDPwAeBrks7tv1NE3BkRzRHR3NRU/fUGEcHqXJ7L5jUxtqG+3OGYWQ3JMllsA2YWbc9I2gayjH63oCJiW/LPLcCDHDueUZPWvbSXV/YeZrHXrjCzEZZlslgDzJM0V1IjhYTwullNkhYCk4DHi9omSRqbvJ8KXAqs779vrWnJ5QtV266vMLMRltlsqIjokXQzcC9QDyyPiHWSbgdaI+Jo4lgG3BURUbT7IuCbkvooJLQvFs+iqlUtbR1cMHMiU121bWYjLLNkARARK4GV/do+22/78wPs9xhwfpaxVZr83sP8on0Pn75mQblDMbMaNFoGuK2Eo2tt+xEfZlYOThYVoqUtz/SJ41l41mnlDsXMapCTRQU43N3LIxu3c+XCaa7aNrOycLKoAI9v3sGh7l5PmTWzsnGyqAAtbR2c0ljPO86ZUu5QzKxGlUwWkk6R9D8l/W2yPU/Sr2cfmsFrVdvveuNUxo1x1baZlUeaXxbfAY4A70y2twH/O7OI7BjrX97LS3sOe+0KMyurNMni3Ij4MtANEBEHAY+yjpDVyUJHly+s/mdfmdnolSZZdEkaDwRA8kC/I5lGZa9a1ZbnrTMnMu20ceUOxcxqWJpk8TngHmCmpO8DLcAfZRqVAdC57wg/37qbJS7EM7MyO+HjPlSY1N8G/BfgHRRuP90aEdtHILaa98DRqm1PmTWzMjthsoiIkLQyIs4H/m2EYrLEqlwHbzhjHOedfXq5QzGzGpfmNtRTkt6eeSR2jKNrbV+5yFXbZlZ+aZ46ewnwQUkvAgco3IqKiHhLppHVuCe27OBgVy+LF3rKrJmVX5pkcU3mUdjrtOTyjB9TzzvPddW2mZVfydtQEfEiMBF4b/KamLRZRiKC1W15LnXVtpmNEmke93Er8H1gWvL6B0m3ZB1YLWt7ZR/bdh9iiWdBmdkokeY21EeASyLiAICkL1FYL/uvsgyslnmhIzMbbdLMhhLQW7Tdix/3kalVuQ7eMuMMpp3uqm0zGx3S/LL4DvCkpJ8m2+8Dvp1dSLVt+/4jPL11N59cPL/coZiZvSrNAPdXgQ8DO5PXhyPia2kOLmmppA2SNkm6bYDP/0LS08nrOUm7iz67QdLG5HVD+lOqbA+05YnACx2Z2ahS8peFpHcA6yLiqWT7dEmXRMSTJfarB+4ArgLagTWSVkTE+qN9IuJTRf1vAS5M3k+m8EyqZgoPMFyb7LtrsCdYaVpyec46fRxveoOrts1s9EgzZvENYH/R9v6krZSLgU0RsSUiuoC7gOtO0P/9wA+T99cA90fEziRB3A8sTfGdFe1ITy8Pb+x01baZjTqpBrgjIo5uREQf6cY6pgNbi7bbk7bXf4E0G5gLrB7MvpJuktQqqbWzszNFSKPbk1t2cqCrl8WeBWVmo0yaZLFF0ickjUletwJbhjmOZcCPI6K3ZM8iEXFnRDRHRHNTU+UvDtSS62DcmDoufePUcodiZnaMNMni48CvUFhOtZ3Cs6JuSrHfNmBm0faMpG0gy3jtFtRg960KEUFLm9faNrPRKc1sqHxELIuIaRFxZkR8ICLyKY69Bpgnaa6kRgoJYUX/TpIWApMoFPoddS9wtaRJkiYBVydtVeu5jv207zrElX5woJmNQmke9/HlZAbUGEktkjolfajUfhHRA9xM4SKfA+6OiHWSbpd0bVHXZcBd/cZFdgJfoJBw1gC3J21Va1WuA/CUWTMbndIMVF8dEX8k6XrgBQqr5j0E/EOpHSNiJbCyX9tn+21//jj7LgeWp4ivKqxuy3P+9DM401XbZjYKpRmzOJpQfg34x4jYk2E8NWnH/iM89ctdfhaUmY1aaX5Z/KukNuAQ8LuSmoDD2YZVWx7c0EkELFnk8QozG53SDHDfRmE2VHNEdAMHOXFxnQ1SS1sH004b66ptMxu10vyyoHhwOXlU+YHMIqoxXT19PPTcdt771rOpq3PVtpmNTmnGLCxD//H8TvYf6fFa22Y2qjlZlNmqXAdjG1y1bWajW5o6i3+S9GuSnFiGWaFqu4NL3ziV8Y2u2jaz0StNAvhr4APARklflLQg45hqxqb8frbuPORCPDMb9dLMhloVER8ELqJQlLdK0mOSPixpTNYBVrNVOa+1bWaVIdWtJUlTgBuBjwL/CXydQvK4P7PIasDqtg7e9IbTOfuM8eUOxczshNKMWfwUeBg4BXhvRFwbET+KiFuAU7MOsFrtOtDF2hd3ee0KM6sIaeos/jIiHhjog4hoHuZ4asYDG/L0BSx21baZVYA0t6HOkzTx6Eby2PD/nmFMNaGlLU/TaWM5f/oZ5Q7FzKykNMniYxGx++hGsib2x7ILqfp19fTx0IZOrlwwzVXbZlYR0iSLekmvXtEk1QON2YVU/da8sJN9R3o8ZdbMKkaaMYt7gB9J+may/d+SNhuillyexoY63jXPVdtmVhnSJIs/ppAgfjfZvh/4VmYRVbmjVdu/cu4UTmlM9RxHM7OyK3m1iog+4BvJy07S5s4DvLjjIB+97Jxyh2JmllrJZCFpHvDnwHnAq2t+RoSvdkPQkqy17aptM6skaQa4v0PhV0UPcAXwd6RYf9sG1pLLs+js05k+0VXbZlY50iSL8RHRAigiXoyIz1NYj7skSUslbZC0SdJtx+nzm5LWS1on6QdF7b2Snk5eK9J832i3+2AXrS/uZIlnQZlZhUkzwnokeTz5Rkk3A9tI8ZiPZIrtHcBVQDuwRtKKiFhf1Gce8Bng0ojYJan4KnooIi4YxLmMeg9u6KQvfAvKzCpPml8Wt1J4LtQngLcBHwJuSLHfxcCmiNgSEV3AXbx+7e6PAXckhX5ERD5t4JVoVa6Dqac28tYZE0t3NjMbRU6YLJJfB78VEfsjoj0iPhwRvxERT6Q49nRga9F2e9JWbD4wX9Kjkp6QtLTos3GSWpP296U5mdGsu7ePnz3XyRWu2jazCnTC21AR0SvpXRl//zzgcmAG8JCk85PHi8yOiG2SzgFWS3omIjYX7yzpJuAmgFmzZmUY5slb88JO9h3u8YMDzawipRmz+M9kgPkfgQNHGyPin0rstw2YWbQ9I2kr1g48GRHdwPOSnqOQPNZExLbke7ZIehC4EDgmWUTEncCdAM3NzZHiXMpmdS5PY30dl7lq28wqUJoxi3HADuBK4L3J69dT7LcGmCdprqRGYBnQf1bTP1P4VYGkqRRuS21Jnmw7tqj9UmA9FaylLc87zp3ChLGu2jazypOmgvvDQzlwRPQks6fuBeqB5RGxTtLtQGtErEg+u1rSeqAX+HRE7JD0K8A3JfVRSGhfLJ5FVWk2d+7n+e0H+PClc8odipnZkKSp4P4O8LpbPBHxO6X2jYiVwMp+bZ8teh/A7yev4j6PAeeXOn6lWO21ts2swqW5J/KvRe/HAdcDL2UTTnValetg4VmnMWPSKeUOxcxsSNLchvpJ8bakHwKPZBZRldlzsJvWF3fx8ff4UVpmVrnSDHD3Nw/w/ZSUHnwuT29fcOVCT5k1s8qVZsxiH8eOWbxCYY0LS6Ell2fKhEYumOmqbTOrXGluQ502EoFUo57ePh7ckOfqN51Fvau2zayClbwNJel6SWcUbU+shsdvjITWF3ex93APiz0LyswqXJoxi89FxJ6jG8mjOD6XXUjVoyXXUajant9U7lDMzE5KmmQxUB+XIafQ0pbnknMmc6qrts2swqVJFq2Svirp3OT1VWBt1oFVuue3H2BL5wHfgjKzqpAmWdwCdAE/orAmxWHg97IMqhocXWvbT5k1s2qQZjbUAWDAJVHt+FpyeeafeSozJ7tq28wqX5rZUPdLmli0PUnSvdmGVdn2HOpmzQs7/avCzKpGmttQU5MZUAAkS6D6RvwJPPRcJz19wZJF/mMys+qQJln0SXp1GTpJsxngKbT2mpZcB5MnNHLBzEnlDsXMbFikmdP5J8Ajkn4GCLiMZClTe72e3j4e2NDJ4kXTXLVtZlUjzQD3PZIuAt6RNH0yIrZnG1bleuqXu9lzqJslHq8wsyqStlqsF8hTWM/iPElExEPZhVW5WnIdjKmX19o2s6qS5qmzHwVuBWYAT1P4hfE4hTW5rZ9VuQ4umTuF08aNKXcoZmbDJs0A963A24EXI+IK4EJg94l3qU0vbD/A5s4DXj7VzKpOmmRxOCIOA0gaGxFtwIJsw6pMLW2FtbY9XmFm1SZNsmhPivL+Gbhf0v8DXkxzcElLJW2QtEnSgFXgkn5T0npJ6yT9oKj9Bkkbk9cNab6v3Fa3dTBv2qnMmuKqbTOrLmlmQ12fvP28pAeAM4B7Su0nqR64A7gKaAfWSFoREeuL+swDPgNcGhG7JE1L2idTeAx6M4WajrXJvrsGdXYjaO/hbp7cspOPXDa33KGYmQ27Qa3BHRE/i4gVEdGVovvFwKaI2JL0vwu4rl+fjwF3HE0CEZFP2q8B7o+Incln9wNLBxPrSHutatu3oMys+gwqWQzSdGBr0XZ70lZsPjBf0qOSnpC0dBD7jiqrc3kmnjKGi2a5atvMqk+5V+VpAOYBl1OYmvuQpPPT7izpJpJq8lmzZpXonZ3evuCBDXmuWOCqbTOrTln+stgGzCzanpG0FWsHVkREd0Q8DzxHIXmk2ZeIuDMimiOiuampfEuXPvXLXew62M1iPzjQzKpUlsliDTBP0lxJjcAyYEW/Pv9M4VcFkqZSuC21BbgXuDp5HPok4OqkbVRqyeVpqBPv9lrbZlalMrsNFRE9km6mcJGvB5ZHxDpJtwOtEbGC15LCegqPFPl0ROwAkPQFCgkH4PaI2JlVrCerJdfBxXMnc7qrts2sSmU6ZhERK4GV/do+W/Q+gN9PXv33XQ4szzK+4fDLHQfZmN/PsovLN2ZiZpa1LG9D1YSWtmStbT/iw8yqmJPFSWrJ5Tm3aQJzpk4odyhmZplxsjgJ+w538+TzO1yIZ2ZVz8niJDy8cTvdveGnzJpZ1XOyOAmrch2cMX4Mb5vtqm0zq25OFkPU2xc8uKGTyxc00VDvP0Yzq26+yg3R01t3sfNAF4s9XmFmNcDJYohacnnq68R7XLVtZjXAyWKIWnJ53j5nEmeMd9W2mVU/J4sh2LrzIBs69nnKrJnVDCeLIVidrLXtKbNmViucLIZgVa6Dc6ZO4JymU8sdipnZiHCyGKT9R3p4cstOr11hZjXFyWKQHtnYSVdvH1cu9HiFmdUOJ4tBWpXLc/q4BprnuGrbzGqHk8Ug9PUFD7TluXzBNMa4atvMaoiveIPwdPtudhzo8niFmdUcJ4tBaMl1uGrbzGqSk8UgtOTyvG32JCae0ljuUMzMRpSTRUrtuw7S9so+lvgWlJnVoEyThaSlkjZI2iTptgE+v1FSp6Snk9dHiz7rLWpfkWWcaRyt2vZTZs2sFjVkdWBJ9cAdwFVAO7BG0oqIWN+v648i4uYBDnEoIi7IKr7BasnlmTPlFM7xWttmVoOy/GVxMbApIrZERBdwF3Bdht+XmQNHenh88w4WLzoTSeUOx8xsxGWZLKYDW4u225O2/n5D0i8k/VjSzKL2cZJaJT0h6X0ZxlnSI5u209Xb5ymzZlazyj3A/S/AnIh4C3A/8L2iz2ZHRDPwAeBrks7tv7Okm5KE0trZ2ZlZkC25Dk4b18Db50zO7DvMzEazLJPFNqD4l8KMpO1VEbEjIo4km98C3lb02bbkn1uAB4EL+39BRNwZEc0R0dzUlE3tQ19fsLqtk/fMb3LVtpnVrCyvfmuAeZLmSmoElgHHzGqSdHbR5rVALmmfJGls8n4qcCnQf2B8RPxi2x627z/iW1BmVtMymw0VET2SbgbuBeqB5RGxTtLtQGtErAA+IelaoAfYCdyY7L4I+KakPgoJ7YsDzKIaES25DuoEl893sjCz2pVZsgCIiJXAyn5tny16/xngMwPs9xhwfpaxpbUql6d59mQmTXDVtpnVLt+EP4GXdh8i9/JervQtKDOrcU4WJ9CSVG37ER9mVuucLE5gda6D2VNO4VyvtW1mNc7J4jgOdvXw6OYdXLlwmqu2zazmOVkcxyMbt9PV08cSPzjQzMzJ4nhWt+U5bayrts3MwMliQH19QUtbnnfPb6KxwX9EZma+Eg7g2Zf20LnPVdtmZkc5WQxgVS5fqNpe4GRhZgZOFgNqyXVw0axJTHbVtpkZ4GTxOq/sOcy6l/Z6+VQzsyJOFv20tHUAeLzCzKyIk0U/Lbk8MyePZ940V22bmR3lZFHkUFcvj27azuKFXmvbzKyYk0WRRzdt50iP19o2M+vPyaJIS1ueCY31XDJ3SrlDMTMbVZwsEhHB6rYOV22bmQ3AV8XEs9v20rH3iKfMmpkNwMki0dLWgQRXLGgqdyhmZqOOk0WiJZfnwpkTmXLq2HKHYmY26mSaLCQtlbRB0iZJtw3w+Y2SOiU9nbw+WvTZDZI2Jq8bsoyzY+9hntm2x7egzMyOoyGrA0uqB+4ArgLagTWSVkTE+n5dfxQRN/fbdzLwOaAZCGBtsu+uLGJdnay17SmzZmYDy/KXxcXApojYEhFdwF3AdSn3vQa4PyJ2JgnifmBpRnHSkutg+sTxLDjztKy+wsysomWZLKYDW4u225O2/n5D0i8k/VjSzEHue9IOd/fyyKbtLFnktbbNzI6n3APc/wLMiYi3UPj18L3B7CzpJkmtklo7OzuHFMDeQ91cfd5ZLH3z2UPa38ysFmSZLLYBM4u2ZyRtr4qIHRFxJNn8FvC2tPsm+98ZEc0R0dzUNLQpr9NOH8dfvv9C3nmuq7bNzI4ny2SxBpgnaa6kRmAZsKK4g6Tiv85fC+SS9/cCV0uaJGkScHXSZmZmZZDZbKiI6JF0M4WLfD2wPCLWSbodaI2IFcAnJF0L9AA7gRuTfXdK+gKFhANwe0TszCpWMzM7MUVEuWMYFs3NzdHa2lruMMzMKoqktRHRXKpfuQe4zcysAjhZmJlZSU4WZmZWkpOFmZmV5GRhZmYlVc1sKEmdwIsncYipwPZhCqecquU8wOcyWlXLuVTLecDJncvsiChZ1Vw1yeJkSWpNM31stKuW8wCfy2hVLedSLecBI3Muvg1lZmYlOVmYmVlJThavubPcAQyTajkP8LmMVtVyLtVyHjAC5+IxCzMzK8m/LMzMrCQni4SkLyQr9j0t6T5Jbyh3TEMl6f9IakvO56eSJpY7pqGS9F8lrZPUJ6niZq5IWippg6RNkm4rdzwnQ9JySXlJz5Y7lpMhaaakByStT/7burXcMQ2VpHGS/kPSz5Nz+V+ZfZdvQxVIOj0i9ibvPwGcFxEfL3NYQyLpamB18pj4LwFExB+XOawhkbQI6AO+CfxhRFTMo4Ul1QPPAVdRWBp4DfD+iFhf1sCGSNK7gf3A30XEm8sdz1Al6+icHRFPSToNWAu8rxL/vaiwFvSEiNgvaQzwCHBrRDwx3N/lXxaJo4kiMQGo2CwaEfdFRE+y+QSFlQYrUkTkImJDueMYoouBTRGxJSK6gLuA68oc05BFxEMU1p2paBHxckQ8lbzfR2HRtenljWpoomB/sjkmeWVy7XKyKCLpzyRtBT4IfLbc8QyT3wH+vdxB1KjpwNai7XYq9KJUrSTNAS4EnixvJEMnqV7S00AeuD8iMjmXmkoWklZJenaA13UAEfEnETET+D5wc3mjPbFS55LydkkxAAACIklEQVT0+RMKqxB+v3yRlpbmXMyGm6RTgZ8An+x3Z6GiRERvRFxA4Q7CxZIyuUWY2bKqo1FELEnZ9fvASuBzGYZzUkqdi6QbgV8HFscoH5gaxL+XSrMNmFm0PSNpszJL7u//BPh+RPxTueMZDhGxW9IDwFJg2Cch1NQvixORNK9o8zqgrVyxnCxJS4E/Aq6NiIPljqeGrQHmSZorqRFYBqwoc0w1LxkU/jaQi4ivljuekyGp6ehsR0njKUymyOTa5dlQCUk/ARZQmHnzIvDxiKjIvwVK2gSMBXYkTU9U8Myu64G/ApqA3cDTEXFNeaNKT9KvAl8D6oHlEfFnZQ5pyCT9ELicwhNOO4DPRcS3yxrUEEh6F/Aw8AyF/98B/kdErCxfVEMj6S3A9yj891UH3B0Rt2fyXU4WZmZWim9DmZlZSU4WZmZWkpOFmZmV5GRhZmYlOVmYmVlJThZmGZJ0lqS7JG2WtFbSSknzyx2X2WDVVAW32UhKir9+CnwvIpYlbW8FzqTwNFqziuFkYZadK4DuiPibow0R8fMyxmM2ZL4NZZadN1NYK8Gs4jlZmJlZSU4WZtlZB7yt3EGYDQcnC7PsrAbGSrrpaIOkt0i6rIwxmQ2Jk4VZRpJ1RK4HliRTZ9cBfw68Ut7IzAbPT501M7OS/MvCzMxKcrIwM7OSnCzMzKwkJwszMyvJycLMzEpysjAzs5KcLMzMrCQnCzMzK+n/A8UyF2aJr3RAAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# scores_：dict with classes as the keys, and the values as the grid of scores obtained during cross-validating each fold,\n",
    "# Each dict value has shape (n_folds, len(Cs))\n",
    "Cs = [1e-3, 1e-2, 1e-1, 1, 10, 100, 1000]\n",
    "n_Cs = len(Cs)\n",
    "#Scores only 1 row to show the \n",
    "scores =  np.zeros((1,n_Cs))\n",
    "scores[0][:] = np.mean(lrcv_L1_roc.scores_[1],axis = 0)\n",
    "    \n",
    "mse_mean = np.mean(scores, axis = 0)\n",
    "plt.plot(np.log10(Cs), mse_mean.reshape(n_Cs,1)) \n",
    "#plt(np.log10(reg.Cs)*np.ones(3), [0.28, 0.29, 0.30])\n",
    "plt.xlabel('C')\n",
    "plt.ylabel('accuracy score')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([0.5       , 0.7908508 , 0.83727952, 0.8392369 , 0.83919986,\n",
       "       0.83908805, 0.83905101])"
      ]
     },
     "execution_count": 27,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "mse_mean"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "1 0.8392368972746331\n"
     ]
    }
   ],
   "source": [
    "best_C = np.argmax(mse_mean)\n",
    "best_score = np.max(mse_mean)\n",
    "print(Cs[best_C], best_score)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 29,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[ 0.40688473,  1.11915174, -0.09953927,  0.03395946, -0.06094794,\n",
       "         0.63694822,  0.2977394 ,  0.12254205,  0.12624947]])"
      ]
     },
     "execution_count": 29,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "lrcv_L1_roc.coef_"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 保存模型，用于后续测试"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 30,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pickle as cPickle\n",
    "\n",
    "cPickle.dump(grid.best_estimator_, open(\"pima_indians_diabetes_L1_org.pkl\", 'wb'))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# SGDClassifier 模型尝试"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 31,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.linear_model import SGDClassifier\n",
    "lr = SGDClassifier(average=True, max_iter=100)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 32,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "accuracy of each fold is:  [0.75974026 0.72727273 0.76623377 0.79738562 0.77124183]\n",
      "cv accuracy is: 0.7643748408454292\n"
     ]
    }
   ],
   "source": [
    "# import warnings filter\n",
    "from warnings import simplefilter\n",
    "# ignore all future warnings\n",
    "simplefilter(action='ignore', category=FutureWarning)\n",
    "\n",
    "# 交叉验证用于评估模型性能和进行参数调优（模型选择）\n",
    "#分类任务中交叉验证缺省是采用StratifiedKFold\n",
    "#数据集比较大，采用5折交叉验证\n",
    "from sklearn.model_selection import cross_val_score\n",
    "accuracy = cross_val_score(lr, X_train, y_train, cv=5, scoring='accuracy')\n",
    "#%timeit loss_sparse = cross_val_score(lr, X_train_sparse, y_train, cv=3, scoring='neg_log_loss')\n",
    "print('accuracy of each fold is: ',accuracy)\n",
    "print('cv accuracy is:', accuracy.mean())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "综合以上测试，数据预处理也会影响数据结果，比如"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 2",
   "language": "python",
   "name": "python2"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.16"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
